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I  ABSTRACT 

In  the  last  year,  the  unsteady,  three-dimensional,  incompressible,  viscous  flow  inter¬ 
actions  between  a  pair  of  vortex  tubes  advected  by  a  uniform  free  stream  and  a  spherical 
particle  held  fixed  in  space  was  investigated  numerically  for  a  range  of  particle  Reynolds 
numbers  between  20  and  100.  Useful  correlations  of  lift  coefficient,  moment  coefficient, 
and  drag  coefficient  with  velocity  fluctuation,  Reynolds  number,  offset  distance,  amd  ini¬ 
tial  vortex  size  have  been  obtained  and  reported.  A  new  mechanism  based  upon  droplet 
lift  has  been  suggested  for  the  dispersion  of  sprays.  Interactions  with  an  array  of  vor¬ 
tices  have  resulted  in  predictive  capabilities  useful  for  the  study  of  droplet  dispersion  in  a 
spray.  The  case  of  large  vortex-to-droplet  diameter  has  been  studied  showing  that  droplet 
deflection  increases  with  vortex  size.  The  investigation  for  the  heat  and  mass  transfer  of 
a  droplet  interacting  with  vortex  tubes  showed  significant  sensitivity  to  the  configuration 
of  the  interaction.  Nusselt  numbers  were  modified  from  their  axisymmetric  values. 

II  OBJECTIVES 

The  objectives  of  this  research  program  axe  to  investigate  the  interactions  of  vaporiz¬ 
ing  droplets  with  a  turbulent  field  of  the  type  encountered  in  gas  turbine  combustors.  It  is 
intended  to  develop  predictive  capability  through  the  use  of  correlations.  There  is  interest 
in  the  wide  frequency  range  including  the  important  and  challenging  high-frequency  end 
of  the  energy  spectrum  where  turbulent  length  scales  are  comparable  to  droplet  size.  The 
full  Navier-Stokes  equations  are  numerically  solved  and  a  simple  mathematical  descrip¬ 
tion  for  the  turbulent  velocity  fluctuation  was  employed.  In  the  mathematical  description, 
turbulent-like  fluctuations  were  simulated  in  a  controlled  way  by  introducing  cylindrical 
vortices  which  have  a  length  scale  of  the  order  of  that  of  the  droplet  and  a  strength 
corresponding  to  a  turbulent  velocity  fluctuation.  From  the  calculations,  instantaneous 
lift,  drag,  and  torque  coefficients,  Nusselt  number,  and  Sherwood  number  are  determined. 
Time-averaged  values  of  these  fluctuating  quantities  are  also  determined.  These  quanti¬ 
ties  will  be  used  to  develop  models  for  droplet  dispersion  and  modifications  of  heating 
and  vaporization  rates  due  to  turbulence. 


Ill  SUMMARY  OF  RESEARCH 


Droplet  —  Turbulent  Interaction  Over  a  Wide  Spectral  Range 


A  theoretical/computational  analysis  of  the  interactions  of  droplets  with  vortical  struc¬ 
tures  is  being  conducted.  Both  the  high-frequency  domain  where  length  scales  of  the 
vortices  are  comparable  to  droplet  size  and  the  low-frequency  domain  where  length  scales 
are  much  larger  than  droplet  size  have  been  analyzed.  This  study  is  leading  to  funda¬ 
mental  insights  about  droplet-turbulent  interactions.  To  simulate  the  three-dimensional, 
unsteady  interactions,  a  Navier-Stokes  solver  was  developed.  The  alternating-direction- 
predictor-corrector  scheme  and  a  pressure  correction  equation  are  employed. 
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In  our  earlier  study  (JFM  1995),  we  examined  the  unsteady,  three-dimensional  inter¬ 
actions  between  an  advected  cylindrical  vortex  tube  and  a  fixed  spherical  particle  whose 
diameter  is  of  the  same  order  of  magnitude  as  the  initial  diameter  of  the  vortex.  The 
particle  Reynolds  number  based  on  the  freestream  velocity  and  the  particle  diameter  was 
in  the  range  20  <  .Re  <  100.  The  initial  size  of  the  cylindrical  vortex  tube  was  in  the 
range  0.25  <  <7  <  4,  where  a  is  the  radius  of  the  vortex  tube  normalized  by  the  particle 
radius.  We  found  that  the  maximum  positive  lift  coefficient  and  the  rms  lift  coefficient 
of  the  sphere  are  linearly  proportional  to  the  circulation  of  the  vortex  tube  at  small  val¬ 
ues  of  a.  However,  at  large  values  of  a,  they  are  linearly  proportional  to  the  mciximum 
fluctuation  velocity  induced  by  the  vortex  tube  but  independent  of  a. 

The  above  work  was  extended  to  study  the  interactions  between  a  pair  of  a^lvected 
vortex  tubes  and  a  stationary  spherical  particle.  A  counter-clockwise  rotating  pair  of 
cylindrical  vortex  tubes  were  initially  located  ten  radii  upstream  from  the  center  of  the 
sphere.  The  lift  and  moment  coefficients  of  the  sphere  interacting  with  a  pair  of  vortex 
tubes  as  a  function  of  time  are  nearly  identical,  respectively,  to  those  of  the  sphere  in¬ 
teracting  with  a  single  vortex  tube  if  the  separation  distance  between  the  tube  centers  is 
less  than  2  y/a  vortex  tube  diameter  for  the  lift  coefficient  and  less  than  y/a  vortex  tube 
diameter  for  the  moment  coefficient,  and  the  total  maximum  fluctuation  velocity  due  to 
the  vortex  tubes  is  used  in  the  case  of  a  pair  of  vortex  tubes. 

Our  results  from  a  spherical  particle  interacting  with  a  single  vortex  and  with  a  pair  of 
vortices  can  be  used  to  estimate  the  trajectory  of  a  moving  spherical  particle  interacting 
with  an  array  of  vortex  tubes.  This  extension  provides  a  better  simulation  of  a  droplet 
in  a  turbulent  flow.  A  spherical  particle  is  injected  into  an  array  with  an  infinite  number 
of  counter-clockwise  rotating  vortices  which  are  located  on  the  negative  X-axis  and  are 
separated  from  one  another  by  the  nondimensional  distance  24.  The  particle  trajectory  as 
a  function  of  time  is  computed  by  solving  a  system  of  two  ordinary  differential  equations 
that  are  the  nondimensional  form  of  the  Newton’s  equation  of  motion  in  the  X  and  Y 
directions.  The  gravity  force  is  neglected  in  this  formulation. 

Figure  1  shows  the  trajectories  of  the  sphere  during  the  dimensionless  time  period 
between  0  and  24  for  initial  particle  Reynolds  numbers  100  and  50  with  particle-to-fluid 
density  ratio  200.  The  initial  vortex  size  is  three  times  larger  than  that  of  the  sphere, 
and  the  initial  offset  distance  of  the  sphere  is  zero.  The  initial  maximum  fluctuation 
velocity  of  the  vortex  is  0.4  normalized  by  the  initial  sphere  velocity.  The  sphere  initially 
moves  upward  due  to  the  upwash  of  the  vortex  and  then  moves  downward  due  to  the 
downwash.  The  maximum  positive  deflection  for  the  case  of  Reo  =  50  is  higher  than  that 
of  Rco  =  100. 

Figure  2  shows  two  trajectories  of  the  sphere  which  are  traced  from  the  initial  injection 
until  particle  Reynolds  number  reaches  one  for  initial  particle  Reynolds  numbers  100  and 
50  with  the  same  initial  parameters  as  in  Figure  1.  A  counter-clockwise  rotating  vortex 
tube  produces  not  only  upwash  downstream  of  itself  and  downwash  upstream  of  itself  but 
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it  also  causes  a  shear  flow  across  the  sphere  when  it  passes  the  sphere.  The  combined  eff'ect 
of  the  downwash  and  the  shear  flow  causes  the  magnitude  of  the  maximum  negative  lift 
to  be  greater  than  the  magnitude  of  the  maximum  positive  lift.  Therefore,  the  average  lift 
coefficient  averaged  over  the  time  span  24  (the  interaction  time  with  the  vortex)  is  small 
but  negative  due  to  the  sheax  flow  effect.  This  small  negative  value  of  the  average  lift 
coefficient  becomes  important  when  the  sphere  interacts  with  an  array  of  many  vortices. 
Thus,  the  sphere  travels  upward  only  for  the  short  initicJ  time  period  and  then  moves 
downwaxd  for  the  most  of  the  time  until  it  stops.  The  final  deflection  ratios  defined  by 
\Y}IXf  \  are  1/27  for  the  case  of  Rco  =  100  and  1/23  for  the  case  of  Rto  —  50.  However, 
the  final  deflection  for  the  case  of  Rco  =  100  is  higher  than  that  of  Rco  =  50,  because 
the  sphere  for  the  case  of  of  Rto  =  100  possesses  higher  initial  momentum  and  it  travels 
farther  than  the  sphere  with  Rcq  =  50. 

So  far,  we  examined  the  trajectory  of  a  spherical  particle  interacting  with  an  array 
of  vortices  whose  sizes  are  in  the  Kolmogorov  length  scale  range  and  are  comparable  to 
the  sphere  size.  In  order  to  understand  completely  a  particle  (or  droplet)  motion  in  a 
turbulent  field,  we  must  investigate  the  motion  of  a  particle  interacting  with  a  vortex 
with  large  (integral)  length  scale  which  is  in  the  opposite  end  of  the  spectrum  to  the 
previous  case.  In  the  following  paragraphs,  we  discuss  the  motion  of  a  spherical  particle 
interacting  with  a  large  vortex. 

Initially,  a  spheric<d  particle  is  injected  into  a  counter-clockwise  rotating  large  vortex 
which  is  located  at  the  origin  of  the  coordinates.  The  upstream  flow  boundary  conditions 
relative  to  the  sphere  are  obtained  by  superimposing  the  sphere  velocity  and  the  induced 
velocity  due  to  the  vortex  tube  at  the  location  of  the  sphere,  and  the  pressure  field  due 
to  the  vortex  tube  at  the  location  of  the  sphere  is  imposed  as  the  pressure  boundary 
condition  at  the  computational  outer  bovmdary  relative  to  the  sphere. 

Figure  3  shows  two  trajectories  of  the  sphere  traced  from  the  initial  location  X  =  300 
to  X  =  —300  for  initial  Reynolds  numbers  100  and  50  with  density  ratio  200.  The  half 
circle  in  Figure  4  represents  the  large  vortex  tube.  The  initial  vortex  size  is  200  times  that 
of  the  sphere,  and  the  initicd  offset  distance  of  the  sphere  is  zero.  The  initial  maximum 
fluctuation  velocity  of  the  vortex  is  0.4  normalized  by  the  initial  sphere  velocity.  The 
sphere  initiedly  moves  upward  due  to  the  upwash  of  the  vortex  and  then  moves  downward 
due  to  the  downwash  of  the  vortex.  The  m£iximum  deflection  for  the  case  of  Rto  =  50 
is  higher  than  that  of  Rto  =  100.  The  maximum  deflection  rates  are  1/6  for  the  case  of 
Reo  =  100  and  1/3.4  for  the  caise  of  Reg  =  50.  These  deflection  rates  are  much  larger 
than  those  from  the  case  of  an  array  of  small  vortices  for  the  same  maximum  fluctuation 
velocity.  Furthermore,  we  know  that  a  vortex  with  the  integral  scale  size  produces  higher 
fluctuation  velocity  than  does  a  small  vortex. 

As  expected,  it  is  found  that  the  large  vortex  causes  larger  deflection  of  the  particle 
than  does  an  array  of  small  vortices.  Also,  the  preliminary  run  for  the  large  vortex  showed 
that  our  formulation  and  computer  code  produce  reasonable  and  believable  solutions.  We 
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are  now  in  the  production  stage  and  will  investigate  the  particle  deflection  as  a  function 
of  initial  Reynolds  number,  vortex  size,  density  ratio,  and  maximum  fluctuation  velocity. 


The  Temperature  Field  and  Vaporization 

There  is  an  ongoing  parciUel  study  focusing  on  enhancement  of  heat  transfer  and 
vaporization  of  the  liquid  droplet  due  to  small-scale  turbulence.  To  gain  confidence  in  the 
code  employed,  we  calculated  droplet  Nusselt  number  for  cm  axisymmetric  flow  and  thus 
successfully  benchmarked  our  code  against  published  literature. 

Following  that,  we  have  since  then  included  the  advecting  cylindrical  vortex  tube  in 
the  simulation  and,  analogous  to  our  parallel  fluid  dynamics  study,  traeed  the  changes  in 
the  temperature  field  due  to  four  parameters:  the  vortex  upstream  offset  distance  (do//), 
flow  field  maximum  velocity  fluctuations  represented  by  the  vortex  initial  maximum  radial 
velocity  {vmax)^  vortex  initial  radius  ((t),  and  the  flow  field  Reynolds  number  (Re). 

Variations  in  the  droplet  Nusselt  number  due  to  the  first  three  parameters  have  been 
determined.  The  influence  of  flow  field  Reynolds  number  on  vortex  interaction  with  the 
droplet  heat  transfer  is  currently  under  study.  Next,  we  will  attempt  to  find  correlations 
between  the  observed  changes  and  the  four  parameters  involved.  It  has  been  observed  that 
the  droplet  temperature  field  could  be  appreciably  affected  due  to  its  collision  with  the 
vortex  tube.  Preliminary  calculations  have  shown  that  under  ‘favorable’  initial  conditions 
of  the  vortex,  the  droplet  Nusselt  number  can  change,  e.g.  due  to  the  parameter  (do/f),  by 
up  to  9%  through  its  interaction  with  the  vortex.  This  is  shown  in  Figure  4.  A  compound 
effect  of  parameters  is  yet  to  be  Studied.  Also,  we  anticipate  that,  when  subsequently 
we  model  the  domain  using  an  initially  stratified  temperature  field  (as  opposed  to  a 
currently  initially  uniform  temperature  field),  more  profound  variations  in  the  droplet 
Nusselt  number  will  be  observed  due  to  the  advecting  vortex. 

Our  investigation  has  also  revealed  that  these  increased  values  of  Nusselt  number  are  a 
consequence  of  the  modification  of  the  near  field  flow  configuration;  the  advecting  vortex 
induces  a  secondary  vortical  perturbation  in  the  sphere  wake  which  changes  orientation 
depending  on  the  location  of  the  parent  vortex.  Future  calculations  will  include  a  strat¬ 
ified  gas;  in  such  a  case,  larger  fluctuations  in  Nusselt  number  should  occur.  Vaporizing 
droplets  will  also  be  considered. 
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Fig.  1.  Two  trajectories  of  the  sphere  during  the 
time  period  between  0  and  24  for  initial 
particle  Reynolds  numbers  100  and  50. 


Fig.  2.  Two  trajectories  of  the  sphere  traced 
from  the  initiad  injection  until  particle 
Reynold  number  reaches  one  for  initial 
pairticle  Reynolds  numbers  100  and  50. 


Fig.  3,  Two  trajectories  of  the  sphere  traced 
from  the  initial  location  to  X  =  —300 
for  initial  particle  Reynolds  numbers 
50  and  100. 


Fig.  4.  Variation  of  the  droplet  Nusselt  number 
versus  the  vortex  initial  offset  distance. 
(This  result  is  for  Vmax  =  0.4,  Re  =  100, 
and  cr  =  1.) 
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Astronautics  and  of  the  American  Physical  Society.  Also,  our  results  have  been  presented 
at  the  AFOSR  Contractors  Meeting.  Professor  Elghobashi  recently  gave  an  invited  lec¬ 
ture  at  an  ASME  Gas-Particle  Flows  meeting.  These  presentations  have  led  to  informal 
interactions  with  government  laboratory  researchers  and  industricd  representatives.  Dis¬ 
cussions  about  the  practical  implications  of  this  research  have  been  held  with  Dr.  Melvin 
Rocquemore  at  Wright  Patterson  Air  Force  Aeronautical  Laboratory  and  with  Dr.  Hukam 
Mongia,  formerly  of  Allison  Gas  Turbine  Division  and  currently  of  General  Electric.  While 
formal  technology  transition  hcis  not  yet  occurred,  it  is  generally  agreed  that  the  research 
results  potentially  can  have  broad  and  lasting  impact. 
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The  unsteady,  three-dimensional,  incompressible,  viscous  flow  interactions  between  a 
vortical  (initially  cylindrical)  structure  advected  by  a  uniform  free  stream  and  a 
spherical  particle  held  flxed  in  space  is  investigated  numerically  for  a  range  of  particle 
Reynolds  numbers  20  <  <  100.  The  counter-clockwise  rotating  vortex  tube  is 

initially  located  ten  sphere  radii  upstream  from  the  sphere  centre.  The  flnite-difference 
computations  yield  the  flow  properties  and  the  temporal  distributions  of  lift,  drag,  and 
moment  coefficients  of  the  sphere.  Initially,  the  lift  force  is  positive  owing  to  the 
upwash  on  the  sphere,  then  becomes  negative  owing  to  the  downwash  as  the  vortex 
tube  passes  the  sphere.  Varying  the  size  of  the  vortex  core  (cr)  shows  that  the  r.m.s.  lift 
coefficient  is  linearly  proportional  to  the  circulation  of  the  vortex  tube  at  small  values 
of  <r.  At  large  values  of  cr,  the  r.m.s.  lift  coefficient  is  linearly  proportional  to  the 
maximum  fluctuation  velocity  («„„*)  induced  by  the  vortex  tube  but  independent  of  <r. 
For  intermediate  values  of  cr,  the  r.m.s.  lift  coefficient  depends  on  both  cr  and  (or 
equivalently  both  cr  and  the  circulation).  We  observe  some  interesting  flow  phenomena 
in  the  near  wake  as  a  function  of  time  owing  to  the  passage  of  the  vortex  tube. 


1.  Introduction 

This  paper  is  concerned  with  the  unsteady,  three-dimensional,  incompressible, 
viscous  flow  interactions  between  a  vortical  structure  (initially  cylindrical)  advected  by 
a  uniform  free  stream  and  a  spherical  solid  particle  which  is  held  flxed  in  space.  This 
flow  is  equivalent  to  that  of  a  spherical  particle  moving  along  a  straight  line  and 
traversing  the  vortical  structure  at  constant  velocity.  The  particle  Reynolds  number 
based  on  the  free-stream  velocity  and  particle  diameter  is  in  the  range  20  <  Re  <  100. 
We  obtain  the  unsteady  velocity  and  pressure  distributions  via  the  numerical  solution 
of  the  time-dependent  three-dimensional  Navier-Stokes  equations  within  a  spherical 
domain  surrounding  the  sphere  and  the  moving  vortex  tube. 

The  motivation  for  studying  this  flow  is  the  need  to  understand  how  the  forces  (drag, 
lift  and  torque)  imparted  on  a  particle  are  influenced  by  fluctuations  in  the  velocity  and 
pressure  of  the  carrier  flow  as  is  the  case  in  particle-laden  turbulent  flows.  Knowledge 
of  the  time-dependence  of  these  forces  and  the  unsteady  flow  field  is  essential  for  the 
accurate  calculation  of  the  particle  trajectory  and  the  heat  and  mass  transfer  rate  of  the 
particle  (or  droplet)  which  in  turn  is  a  prerequisite  for  predicting  particle  dispersion 
and  vaporization  rate  in  turbulent  flows.  The  exact  relations  between  these  forces  and 
the  turbulent  fluctuations  cannot  be  obtained  analytically  owing  to  the  nonlinearity  of 
the  equations  governing  the  motion  of  the  particle  and  fluid.  While  the  cylindrical 
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vortex  is  far  too  simple  to  represent  real  turbulence,  some  important  elementary 
understanding  can  result  from  this  study. 

Numerical  simulation  of  the  dispersion  of  particles  in  a  turbulent  flow  requires  the 
solution  of  the  equation  of  particle  motion.  This  equation  which  is  classically  known 
as  the  Basset-Boussinesq-Oseen  (BBO)  equation  and  has  been  re-derived  recently  by 
Maxey  &  Riley  (1983)  is  restricted  to  low  Reynolds  number  R^^l,  where 
Re  =  d  \u'  —  v'\lv'\  v'  and  «'  are  the  velocities  of  the  particle  and  its  surrounding  fluid 
respectively,  d'  is  the  particle  diameter,  and  i>'  is  the  fluid  kinematic  viscosity. 
Furthermore,  the  drag  force  in  that  equation  consists  of  two  terms,  namely,  the  quasi¬ 
steady  Stokes  drag  and  the  unsteady  memory  term  (Basset).  The  former  is  purely 
viscous,  whereas  the  latter  depends  on  both  the  viscosity  and  particle  acceleration 
relative  to  the  fluid.  The  superposition  of  these  two  terms  is  a  result  of  the  linearization 
of  the  Navier-Stokes  equations  by  Basset  (1888).  A  more  serious  restriction  (than 
Re  -4  1)  in  the  equation  of  particle  motion  is  that  the  velocity  gradients  in  the  carrier 
flow  in  the  neighbourhood  of  the  particle  should  be  very  small.  This  requires  that  the 
shear  Reynolds  number  =  {a'^/v'){U'JL')  1,  where  a'  is  the  particle  radius,  and 

{U'JU)  is  a  reference  gradient  of  the  undisturbed  velocity  field.  Therefore,  the 
interesting  case  in  which  the  eddy  (or  vortex)  size  is  comparable  to  that  of  the  particle 
cannot  be  properly  treated  by  the  standard  equation  of  particle  motion.  This  situation, 
in  addition  to  being  relevant  to  the  fundamental  understanding  of  fluid  dynamics,  is  of 
practical  interest  as  well.  For  example,  in  a  typical  gas  turbine  combustor  where  the 
Reynolds  number  is  of  the  order  of  10*  and  the  integral  lengthscale  is  of  the  order  of 
0.1  m,  the  smallest  (Kolmogorov)  lengthscale,  rj,  is  about  100  p,  which  is  comparable 
to  the  size  of  a  typical  fuel  droplet.  Fluid  motion  at  the  Kolmogorov  lengthscale 
experiences  the  largest  strain  rates  and  scalar  gradients  in  the  flow.  The  largest  scalar 
gradients  control  the  important  phenomena  of  heat  and  mass  transfer  and  chemical 
reaction.  Motion  at  the  largest  lengthscales  ( >  17)  contains  most  of  the  turbulence 
energy  and  governs  the  dispersion  of  particles  (or  droplets)  but  not  the  small-scale 
phenomena  mentioned  above. 

Almost  all  application-oriented  studies  of  diluted  particle  suspension  calculate  the 
drag  on  the  particle  using  the  standard  drag  curve.  This  drag  curve  has  been  obtained 
(experimentally  and  numerically)  for  a  particle  fixed  in  space  subjected  to  a  steady  flow. 
In  the  case  of  unsteady  flow,  this  drag  relationship  is  an  approximation  that  can  be 
valid  only  if  the  timescale  of  the  particle  motion  is  much  larger  than  that  of  the  flow. 
Note  that  empirical  relationships  for  unsteady  drag  have  been  proposed  (e.g. 
Houghton  1963;  Odar  1966;  Schonebom  1975;  Ingebo  1956).  The  expression  of 
Ingebo,  derived  experimentally,  is  valid  only  for  the  limited  conditions  of  the 
experiment.  All  other  expressions  concern  mainly  purely  harmonic  flows.  Additional 
interactions  between  the  particle  and  the  flow  are  the  well-known  Saffman’s  lift  due  to 
uniform  shear  (1965,  1968)  and  the  lift  due  to  particle  rotation  (Rubinow  &  Keller 
1961).  Saffman’s  lift  force  expression  is  valid  only  for  Re  <  R\'hear  ^'*‘1  ^shear  ^ 
where  =  {6u' ldy')d^l{Av').  Under  these  conditions,  the  lift  due  to  particle 

rotation  is  negligible  (Saffman  1965).  Recently,  McLaughlin  (1991)  removed  the 
restriction  Re  R\hear  provided  a  new  form  for  the  lift  force. 

Three-dimensional  flow  interactions  between  a  vortical  structure  and  a  bluff  body  (a 
rigid  sphere  in  the  simplest  form)  at  finite  Reynolds  number  have  not  been  investigated 
yet.  Our  present  approach,  outlined  at  the  beginning  of  this  section,  is  a  first  step 
toward  better  understanding  of  the  physics  of  interaction  between  a  particle  and  the 
carrier  turbulent  flow.  For  example,  we  examine  the  details  of  the  temporal  behaviour 
of  the  flow  structure  around  the  sphere  owing  to  the  passage  of  the  vortex  tube. 
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Figure  1.  Flow  geometry  and  coordinates. 


Furthermore,  we  study  the  effect  of  varying  the  ratio  of  vortex  tube  size  to  particle  size, 
Reynolds  number,  and  offset  distance  between  the  particle  and  the  vortex  tube  on  the 
temporal  distributions  of  the  forces  imparted  on  the  particle  (drag,  lift  and  torque)  and 
the  flow  structure  in  the  neighbourhood  of  the  particle. 

The  next  section  provides  a  mathematical  description  of  the  flow  considered,  the 
governing  equations  and  the  numerical  solution  procedure.  Section  3  discusses  the 
results  including  the  numerical  accuracy  issues  and  the  effects  of  varying  the 
parameters  listed  above.  Section  4  provides  a  summary  of  the  work. 


2.  Problem  statement  and  formulation 

2.1.  The  flow  description 

We  consider  the  time-dependent,  three-dimensional,  incompressible,  viscous  flow 
interactions  between  an  initially  cylindrical  vortex  tube  and  a  spherical  solid  particle. 
The  vortex  tube  is  moving  with  the  laminar  free  stream,  and  a  sphere  is  suddenly  placed 
and  held  fixed  in  space  as  shown  in  figure  1 .  The  initial  offset  distance,  d^ff,  denotes  the 
shortest  distance,  normalized  by  the  sphere  radius,  between  the  initial  vortical  axis  and 
the  (>’,r)-plane,  which  is  parallel  to  the  free  stream.  All  the  variables  are  non- 
dimensionalized  using  the  sphere  radius  a'^  as  the  characteristic  length  and  U'^  as  the 
characteristic  velocity,  where  the  superscript  '  denotes  dimensional  quantity.  The 
cylindrical  vortex  tube,  whose  diameter  is  of  the  order  of  the  sphere  diameter,  is 
initially  located  ten  radii  upstream  from  the  centre  of  the  sphere.  The  effects  of  the 
vortex  tube  on  the  sphere  are  negligible  at  this  initial  distance  because  the  magnitude 
of  the  initial  velocity  field  induced  by  the  vortex  tube  is  less  than  2%  of  the  free-stream 
velocity.  Far  upstream,  the  flow  is  uniform  with  constant  velocity  (7^  Ac  parallel  to  the 
(y,z)-plane.  We  have  one  symmetry  plane,  the  (x,  z)-plane,  as  seen  in  figure  1. 

Two  coordinate  systems  are  used  in  our  formulation:  the  Cartesian  coordinates 
{x,y,z)  and  the  non-orthogonal  generalized  coordinates  0-  The  origin  of  the 
former  coincides  with  the  sphere  centre.  ^  is  the  radial,  7f  is  the  angular,  and  ^  is  the 
azimuthal  coordinates.  The  non-orthogonal  generalized  coordinate  system  can  be 
easily  adapted  to  three-dimensional  arbitrary  geometries.  In  the  present  study,  a 
spherical  domain  is  used,  and  the  grid  reduces  to  an  orthogonal,  spherical  grid.  The 
grids  are  denser  near  the  surface  of  the  spherical  particle,  and  the  grid  density  in  the 
radial  direction  is  controlled  by  the  stretching  function  developed  by  Vinokur  (1983). 
Owing  to  symmetry,  the  physical  domain  is  reduced  to  a  half  spherical  space.  The 
domain  of  the  flow  is  bounded  by  1  <  ^  <  Wj,  1  ^  i;  <  ATj,  1  <  ^  <  V3,  where  ^  =  1 
and  Vj  correspond,  respectively,  to  the  sphere  surface  and  the  far-field  boundary 
surrounding  the  sphere;  rj  =  1  and  denote,  respectively,  the  positive  z-axis  and  the 
negative  z-axis;  ^  =  1  and  refer,  respectively,  to  the  (x,z)-plane  in  the  positive 
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Figure  2.  Comparison  of  tangential  velocities  induced  by - ,  a  point  vortex  and 

- ,  a  vortex  tube  for  r„  =  2.5  and  <r  =  1. 

x-direction  and  the  (x,z)-plane  in  the  negative  x-direction.  Uniform  spacing 
{8i  =  Si}  =  S^=  I)  is  used,  for  convenience,  for  the  generalized  coordinates. 

2.2.  The  vortex  tube  features 

The  initial  vortex  tube  has  a  small  core  region  with  a  radius  o’  (normalized  by  the 
sphere  radius).  This  core  is  defined  such  that  the  initial  velocity  induced  by  the  vortex 
tube  approaches  zero  as  the  distance  from  the  centre  of  the  vortex  tube  goes  to  zero, 
and  at  distances  much  greater  than  o-,  the  induced  velocity  approaches  that  of  a  point 
vortex  (figure  2).  We  use  the  vortex  tube  construction  of  Spalart  (1 982),  which  has  the 
following  stream  function: 


U^,z.t  =  0)  =  -^ln[(.v-.v„)^  +  (r-z„)^  +  ^U 


(I) 


where  T,,  is  the  non-dimensional  circulation  around  the  vortex  tube  at  radius  cr  and  at 
the  initial  time.  is  positive  when  the  vortex  tube  rotates  counterclockwise,  and 
and  Zo  denote  the  location  of  the  centre  of  the  vortex  tube.  The  circulation  around  a 
circular  path  far  away  from  the  centre  of  the  vortex  is  given  by  T,  =  2r„.  The  tangential 
velocity  distribution  of  the  vortex  tube  compared  with  a  point  vortex  is  shown  in  figure 
2  for  Tp  =  2.5  and  a- —  l.O.  As  shown  in  figure  2,  the  cylindrical  vortex  tube  can  be 
viewed  as  an  evolution  from  the  point  vortex  owing  to  the  cylindrical  viscous  diffusion. 
The  initial  pressure  field  due  to  the  vortex  tube  is  obtained  by  solving  the  radial 
component  of  the  Navier-Stokes  equations  which  balances  the  centrifugal  acceleration 
and  the  pressure  gradient  for  circular  streamlines,  and  has  the  following  form 


p^{x,z,t  =  0)  = 


n _ I _ 

2n^  (x  -Xof  +  {z-  zff  -H  o-^  ’ 


(2) 


where  is  non-dimensional  pressure  defined  by  =  ip'„—p'^)/p'U''K  The  pressure 
due  to  the  vortex  tube  attains  its  lowest  value,  p^  =  -  r^/(2nV^),  at  the  centre  of 
the  vortex  tube  and  approaches  zero  at  far  distances  from  the  centre  of  the  vortex. 
Equation  (2)  is  used  to  prescribe  only  the  initial  pressure  field  generated  by  the  vortex 
tube. 
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Figure  3.  (a)  Velocity  and  (b)  vorticity  fields  due  to  a  vortex  tube  as  a  function  of  radial  distance 
and  time  for  Re  =  100,  =  2.5  and  <r  =  1. 


In  order  to  gain  insight  about  the  properties  of  the  vortex  tube,  we  examine  the  flow 
field  generated  in  the  absence  of  the  particle.  We  compute  the  induced  velocity  and 
vorticity  field  as  a  function  of  radius  and  time  owing  to  the  vortex  tube  moving  with 
the  free  stream  (U'^k).  The  origins  of  the  moving  coordinate  system  is  the  centre  of  the 
vortex  tube.  We  solve  the  following  linear  diffusion  equation  which  is  the  tangential 
component  of  the  Navier-Stokes  equations  balancing  the  unsteady  and  diffusion  terms 
for  the  tangential  momentum  (Batchelor  1967), 


dt 


2  1 

1  dUff 

Re' 

RdR 

R'\ 

(3) 


where  R  is  the  radial  distance  from  the  centre  of  the  vortex  tube,  is  the  tangential 
velocity  around  the  vortex  tube  normalized  by  the  free-stream  velocity,  and  Re  is  the 
Reynolds  number  based  on  the  reference  lengthscale  a'„  and  the  free-stream  velocity. 
Figures  3  (a)  and  3  (b)  show,  respectively,  the  velocity  and  vorticity  fields  as  a  function 


128  /•  Kim,  S.  Elghobashi  and  W.  A.  Sirignano 

of  radial  distance  and  time  for  Re  =  100.  F,.  =  2.5,  and  <r  =  1.0.  The  size  of  the  vortex 
core  becomes  larger  as  time  increases  owing  to  viscous  diffusion,  whereas  the 
magnitudes  of  the  tangential  velocity  and  the  vorticity  inside  the  vortex  core  decrease. 
Note  that  this  classical  linearized  analysis  is  not  employed  in  the  present  study;  rather 
a  fully  nonlinear  computational  analysis  is  performed. 

2.3.  Governing  equations  and  boundary  conditions 

The  continuity  and  momentum  equations  and  the  initial  and  boundary  conditions  are 
non-dimensionalized  using  the  sphere  radius  a#  as  the  characteristic  length  and  V'^  as 
the  characteristic  velocity. 

V-K=0.  (4) 

^+V.(Fn  =  -Vp+J-V*K  (5) 

ot  Re 


The  governing  equations  (4)  and  (5)  are  cast  in  terms  of  the  generalized  coordinates 
ii,  1?,  0  lo  a  three-dimensional  body  of  arbitrary  shape.  The  numerical  integration 
is  performed  using  a  cubic  computational  mesh  with  equal  spacing  =  Hrj  =  S^=  1). 

The  velocities  on  the  sphere  surface  are  zero  owing  to  the  no-slip  condition,  and  the 
pressure  on  the  sphere  is  obtained  from  the  momentum  equation.  The  boundary 
conditions  are 


0«  Re  0«* 


M  =  i;  =  W  =  0 


at^=  1, 


(6) 


p  =  0,  u  =  v  =  0,  w=l  at  ^  ^2  (upstream),  (7) 


p  =  0,  ^  ^  =  0  at  ^  \  (downstream),  (8) 


0U 


u  =  0  at  ^  =  1, 


(9) 


where  «,  v,  and  w  are  the  velocities  in  the  .v-,  y-,  and  z-directions,  respectively,  V„  is  the 
velocity  in  the  direction  normal  to  the  sphere  surface,  and  p  is  the  pressure,  n  denotes 
the  direction  normal  to  the  sphere  surface,  d/dn  =  (^* -l-|J-l-f|)*'*0/0|,  and  q  = 
denotes  the  mid-plane  between  rj  =  1  and  iV,.  Equation  (9)  corresponds  to  the 
symmetry  conditions  and  zero  v  velocity  in  the  (x,z)  symmetry  plane.  Conditions 
guaranteeing  continuity  in  the  rj  direction  for  p,  u,  v  and  w  on  the  axes  p  =  1  and  q  = 
are  also  imposed. 

In  order  to  start  the  numerical  solution  of  equations  (4)  and  (5),  we  provide  initial 
velocity  and  pressure  fields  by  superposing  the  flow  fields  due  to  the  uniform  stream 
and  the  vortex  tube  in  addition  to  the  no-slip  condition  on  the  sphere  surface: 


Po=Pv<  except  at  ^=1,  (10) 

Po  =  Pv^  Mo  =  =  ‘*’0  =  0  at  ^  =  1,  (1 1) 

where  and  p^  are  given  by  equations  (1)  and  (2),  respectively. 

The  only  non-dimensional  groupings  appearing  in  the  equations  and  initial  and 
boundary  constraints  are  the  sphere  Reynolds  number,  vortex  tube  radius,  offset 
distance,  and  vortex  circulation  (or  vortex  Reynolds  number). 
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The  drag,  lift  and  moment  coefficients  are  evaluated  in  dimensional  form  as  follows. 


— + j  'kdS\ 

(12) 

—p'n*idS'  +  n*T*idS\ 

(13) 

s  Js 

M'  =  j  r'xt'dS', 

(14) 

where  S'  denotes  the  surface  of  the  sphere,  n  is  the  outward  unit  normal  vector  at  the 
surface,  r'  is  the  position  vector  from  the  centre  of  the  sphere,  and  t'  is  the  viscous  stress 
tensor.  The  lift  force  is  assumed  positive  when  it  is  directed  toward  the  positive  x-axis. 
Owing  to  symmetry,  only  the  ,v-component  of  the  moment  is  non-zero  and  is  assumed 
positive  in  counter-clockwise  direction. 

The  non-dimensional  coefficients  of  drag,  lift  and  moment  are  defined,  respectively, 
as 


f' 

^’>-\p'U'^nar 

(15) 

C^  =  ^p'U'^na',\ 

2 

(16) 

^  Af'.y 

\p'U'iTia,^- 

(17) 

Note  that  in  this  analysis,  the  sphere  does  not  accelerate  or  rotate  owing  to  the 
aerodynamic  forces  and  torque. 


2.4.  Numerical  solution 

We  have  developed  a  three-dimensional,  implicit,  finite-difference  algorithm  to  solve 
simultaneously  the  set  of  the  discretized  partial  differential  equations.  The  method  is 
based  on  an  Altemating-Direction-Predictor-Corrector  (ADPC)  scheme  to  solve  the 
time-dependent  Navier-Stokes  equations.  ADPC  is  a  slight  variation  of  Altemating- 
Direction-Implicit  (ADI)  method.  It  is  first-order  accurate  in  time  but  is  effective  and 
implemented  easily  when  embedded  in  a  large  iteration  scheme  (Patnaik  1986;  Patnaik 
et  al.  1986).  The  control  volume  formulation  is  used  to  develop  the  finite-difference 
equations  from  the  governing  equations  with  respect  to  the  generalized  coordinates 
a,  1/,  Q.  One  of  the  advantages  of  the  control  volume  formulation  is  that  mass  and 
momentum  are  conserved  over  a  single  control  volume,  and  hence  the  whole  domain 
regardless  of  the  grid  fineness.  An  important  part  of  solving  the  Navier-Stokes 
equations  in  primitive  variables  is  the  calculation  of  the  pressure  field.  In  the  present 
work,  a  pressure  correction  equation  is  employed  to  satisfy  indirectly  the  continuity 
equation  (Anderson,  Tannehill  &  Pletcher  1984).  The  pressure  correction  equation  is 
of  the  Poisson  type  and  is  solved  by  the  Successive-Over-Relaxation  (SOR)  method. 

The  overall  solution  procedure  is  based  on  a  cyclic  series  of  guess-and-correct 
operations.  The  velocity  components  are  first  calculated  from  the  momentum  equations 
using  the  ADPC  method,  where  the  pressure  field  at  the  previous  timestep  is  employed. 
This  estimate  improves  as  the  overall  iteration  continues.  The  pressure  correction  is 
calculated  from  the  pressure  correction  equation  using  the  SOR  method,  and  new 
estimates  for  pressure  and  velocities  are  obtained.  This  process  continues  until  the 
solution  converges  at  each  time-step. 


130 


/.  Kim,  S.  Elghobashi  and  W.  A.  Sirignano 


C  C 

Re  =  20 

Cd 

C* 

21x21  x2l 

1.048 

\J11 

2.825 

31  x31  x31 

1.039 

1.740 

2.779 

41  x4I  x41 

1.037  1.731 

Re  =  100 

2.768  , 

2.74 

21x21x21 

0.560 

0.590 

1.150 

31  x31  x3i 

0.535 

0.582 

1.117 

41  x41  x41 

0.526 

0.581 

1.107 

1.09 

Table  1.  Drag  coefficients  as  a  function  of  grid  density  at  Re  =  20  and  100,  where  ♦  denotes  the 
data  from  the  correlation  of  Clift  el  al.  (1978). 


3.  Results  and  discussion 

In  §(3.1),  we  test  the  accuracy  of  the  full  three-dimensional  solution  procedure  by 
predicting  the  axisymmetric  flow  over  a  single  sphere  and  by  examining  the  effects  of 
grid  resolution  on  the  maximum  lift  coefficient  of  the  sphere  owing  to  the  interaction 
between  a  vortex  tube  and  a  sphere.  In  §§(3.2),  (3.3),  (3.4)  and  (3.5),  we  discuss  the  three- 
dimensional  interactions  between  a  vortex  tube  and  a  sphere,  the  effects  of  the  offset 
distance,  the  size  of  the  vortex  tube,  and  Reynolds  number,  respectively. 

3.1.  Numerical  accuracy 

Here  we  examine  the  flow  generated  by  an  impulsively  started  solid  sphere  in  a 
quiescent  fluid  at  two  Reynolds  numbers:  20  and  100.  The  time-dependent  solution 
converges  asymptotically  to  a  steady-state  which  is  in  good  agreement  with  the 
available  experimental  data  and  correlations  as  shown  in  tables  1  and  2.  Table  1  lists 
the  drag  coefficients  as  a  function  of  the  computational  grid  density  at  Reynolds 
numbers  20  and  100,  respectively,  and  compares  them  with  the  correlations  of  Clift, 
Grace  and  Weber  (1978).  Table  2  shows  the  pressures  at  the  front  and  rear  stagnation 
points  and  the  separation  angle  measured  from  the  front  stagnation  point  as  a  function 
of  grid  density  at  Reynolds  number  20  and  100,  in  comparison  with  the  data  of  Taneda 
(1956)  and  also  with  the  correlations  of  Clift  et  al.  (1978).  Although  the  solutions  in 
these  test  cases  are  axisymmetric,  none  of  the  three  velocity  components  in  our 
formulation  becomes  identically  zero.  Therefore,  the  three-dimensional  solution 
scheme  is  fully  exercised  here.  The  calculations  were  performed  for  three  different 
grids,  =  (31x31x31)  and  (41x41x41),  in  a  comp¬ 

utational  domain  with  an  outer  boundary  located  at  21  sphere  radii  from  the  sphere 
centre. 

We  tested  the  solution  procedure  by  varying  the  far-lield  boundary  condition  and  by 
changing  the  location  of  the  outer  boundary.  In  the  first  test,  the  far-field  outflow 
boundary  condition  was  changed  from  0^/0r  =  0  to  0^/0z  =  0  (^  =  u,v  or  w).  There 
was  almost  no  difference  in  the  drag  coefficient  and  the  near  wake  size  (the  separation 
angle  and  length  of  the  recirculation  eddy)  at  Reynolds  numbers  20  and  100.  Our 
calculation  shows  that  separation  does  not  occur  at  Reynolds  number  20.  In  the  second 
test,  the  location  of  the  outer  boundary  in  downstream  was  changed  from  21  to  41 
sphere  radii.  There  was  virtually  no  change  in  the  drag  coefficient  and  the  near  wake 
size  at  both  Reynolds  numbers. 

We  examined  the  effects  of  grid  resolution  on  the  lift  coefficients  of  the  sphere  owing 
to  the  flow  interaction  between  a  cylindrical  vortex  tube  flowing  with  the  free  stream 
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Figure  4.  Pressure  and  shear  stress  distributions  around  the  sphere  in  the  (.r,  z)-plane  of  symmetry 
as  a  function  of  grid  resolution  at  r  =  \2  for  Re  =  100,  =  0,  rr  =  1,  and  =  0.4. 


yVj  X  ^  jVg 

^or 

11 

21  x21  x21 

0.668 

-0.201 

180 

31x31x31 

0.707 

-0.178 

180 

41  x41  x41 

0.721 

-0.171 

Re  =  100 

180 

180 

21  x2l  x21 

0.554 

-0.093  8 

124.09 

31  x31  x31 

0.554 

-0,0831 

125.62 

41  x41  x41 

0.554 

-0.0792 

126.16 

126.5 

Table  2.  Pressure  at  the  front  and  rear  stagnation  points  and  the  separation  angle  measured  from  the 
front  stagnation  point  as  a  function  of  grid  density  at  Re  =  20  and  100,  where  *  denotes  the  data  from 
Taneda  (1956)  and  the  correlation  of  Clift  et  al.  (1978). 


iVj  X  A^2  X  N.^ 

C  C 

^L.max'l.v 

Re  =  20 

max  2 

21  x21  x2I 

-0.336 

-0.359 

-0.695 

31  x31  x31 

-0.358 

-0.363 

-0.722 

41  x41  x41 

-0.362  -0.363 

Re  =  100 

-0.725 

21  x21  x21 

-0.336 

-0.159 

-0.495 

31  x31  x31 

-0.433 

-0.186 

-0.619 

41  x4I  x41 

-0.476 

-0.192 

-0.668 

51  x5I  x51 

-0.494 

-0.191 

-0.685 

Table  3.  Maximum  negative  lift  coefficients  as  a  function  of  grid  density  at  Re  =  20  and  100. 


and  a  sphere  fixed  in  space  at  Reynolds  numbers  20  and  100.  The  lift  coefficients  are 
obtained  for  the  offset  distance  =  vortex  core  radius  cr=  1,  and  maximum 
fluctuation  velocity  =  0.4.  Table  3  shows  the  maximum  negative  lift 

coefficient  (Q;„ax2)  ^he  sphere  as  a  function  of  the  computational  grid  density  at 
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Figure  5.  Pseudo-streamlines  (left-hand  column)  and  contour  lines  of  >’-component  vorticity  (right- 
hand  column)  in  the  principal  plane  at  (a)  t  =  0,  (6)  1,  (c)  6,  (d)  9,  (e)  10,  (f)  1 1,  (g)  12,  (h)  13,  (0 
15,  (;)  18,  (k)  21  and  (/)  30  for  Re  =  100,  =  0,  tr  =  1  and  t)„„  =  0.4. 
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Reynolds  numbers  20  and  100.  The  calculations  were  performed  for' three  different 
grids.  (A^,  X  N,y.N^)  =  (21  x  21  x21),  (31x31x31)  and  (41  x  41  x  41)  for  Re  =  20, 
and  four  different  grids.  (iV,  xN^y.  N^)  =  (21  x 21  x 21),  (31  x  31  x 31),  (41  x 41  x 41) 
and  (51  X  51  X  51)  for  Re  =  100,  in  a  computational  domain  with  an  outer  boundary 
located  at  21  sphere  radii  from  the  sphere  centre.  The  result  of  the  31  x  31  x  31  grid 
differs  by  0.42%  for  that  of  the  41  x41  x41  grid  for  Re  =  20,  and  the  result  of  the 
41  x41  x41  grid  differs  by  2.5  %  from  that  of  the  51  x  51  x  51  grid  for  Re  =  100.  Figure 
4  provides  additional  results  on  the  effect  of  grid  resolution  on  convergence  and  shows 
the  distributions  of  the  pressure  and  shear  stress  coefficients  (normalized  by  the 
dynamic  pressure)  around  the  sphere  in  the  (.v,  z)-plane  of  symmetry  in  the  positive  x- 
direction  for  the  same  parameters  as  used  above  with  Re  =  100.  The  pressure  and  shear 
stress  distributions  were  obtained  at  /  =  12  about  which  the  lift  coefficient  reaches  its 
maximum  in  negative  value.  The  pressure  coefficient  at  (x,z)  =  (0,  —  1)  of  the 
41x41x41  grid  differs  by  0.93%  from  that  of  the  51x51x51  grid.  The  same 
calculations  were  performed  by  changing  the  location  of  the  outer  boundary  in 
downstream  from  21  to  41  sphere  radii.  There  was  virtually  no  change  in  the  lift, 
moment  and  drag  coefficients. 

In  order  to  examine  the  far-field  boundary  effects,  we  repeated  the  simulation  as 
above  for  Re  =  100,  d„ff  =  0,  <r  =  1,  and  v„ax  =  0.4  but  with  a  box-type  computational 
domain  with  symmetry  boundary  conditions  on  its  sides.  The  lift,  moment  and  drag 
coefficients  of  the  box-type  computational  domain  at  /  =  12  differs  by  0.12%,  0.17% 
and  0.13  %,  respectively,  from  those  of  the  spherical  computational  domain  used  in  the 
present  paper.  The  spherical  computational  domain  gives  a  little  finer  resolution  than 
does  the  box-type  computational  domain  with  the  same  number  of  grid  points,  and  so 
smoother  contour  lines  for  the  vorticity  and  stream  lines  in  the  {x,y)  symmetry  plane. 
We  also  solved  the  same  problem  as  above  by  employing  a  complete  computational 
domain  without  the  symmetry  plane  and  periodic  boundary  condition  in  the  ^- 
direction.  In  that  case,  the  lift  coefficients  differ  by  0.16%  from  those  in  table  3  where 
the  symmetry  condition  was  employed.  The  41x41x41  grid  is  used  in  the  following 
calculations. 

The  run  for  the  interaction  between  a  vortex  tube  and  a  sphere  at  Reynolds  number 
100  with  the  41  x41  x41  grid  required  2.62  mega  words,  a  dimensionless  timestep  of 
A/  =  0.0025,  and  a  total  time  of  2.95  c.p.u.  hours  on  Cray  Y-MP8/864  for  the  final 
time  of  tf  =  lA.  Each  timestep  takes  about  1.11  c.p.u.  seconds.  Another  test  was 
performed  to  examine  the  effects  of  the  timestep.  The  same  calculation  as  above  was 
repeated  for  Re  -  100,  d^ff  =  0,  cr  =  1,  and  =  0.4  but  with  the  time-step  reduced 
by  half.  The  lift,  moment  and  drag  coefficients  differ  by  0.25%,  0.19%  and  0.003%, 
respectively,  from  those  obtained  with  the  timestep  used  in  the  presented  results. 

3.2.  Interactions  of  a  vortex  tube  and  a  sphere 

We  consider  the  interactions  of  a  vortical  structure  advected  by  the  free  stream  and  a 
sphere  suddenly  placed  in  the  flow  and  held  fixed  in  space  at  Reynolds  number  100.  The 
vortical  structure  is  initially  a  cylindrical  vortex  tube  rotating  counter-clockwise  in 
figure  1  with  a  non-dimensional  radius  of  unity  and  an  offset  distance  of  zero,  and 
located  at  10  sphere-radii  upstream  from  the  centre  of  the  sphere. 

3.2.1.  Flow  structure 

In  order  to  illustrate  better  the  fluid  motion,  we  consider  the  flow  field  in  the  (x,  z)- 
plane  of  synunetry,  which  is  defined  as  the  principal  plane,  where  the  strongest 
interactions  occur  between  the  vortical  structure  and  the  sphere. 


(0  ij) 


Angle  (deg.)  Angle  (deg.) 


Figure  6.  Pressure  and  shear  stress  distributions  around  the  sphere  in  the  principal  plane  at  (a)  /  =  1, 
(b)  6,  (c)  9,  (d)  10,  (e)  11,  (/)  12,  (g)  13,  (h)  15,  (/)  18  and  {J)  21  for  Re  =  100,  =  0,  <7  =  1  and 
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Pseudo-streamlines  are  employed  in  the  following  illustrations.’  The  pseudo¬ 
streamlines  are  obtained  from  the  pseudo-stream  function  which  is  defined  by 
assuming  that  the  velocity  field  in  the  principal  plane  does  not  change  in  the 
perpendicular  direction  to  the  principal  plane  and  by  using  the  two-dimensional  stream 
function  definition.  The  sphere  surface  in  the  principal  plane  is  used  as  a  reference 
streamline  =  0).  We  note  that  a  real  stream  function  \jr  cannot  be  defined  and 
calculated  from  the  velocity  in  the  principal  plane  owing  to  the  existence  of  a 
divergence  associated  with  the  third  component  of  velocity.  Nevertheless,  for 
descriptive  purposes  only,  it  is  convenient  to  use  the  two-dimensional  stream  function 
definition  to  present  approximations  to  the  streamline  pattern. 

Figures  5  (a)-5(/)  display  the  pseudo-streamlines  (left-hand  column)  and  the  contour 
lines  of  y-component  vorticity  (right-hand  column)  in  the  principal  plane  at  /  =  0,  1 , 
6,  9,  10,  11,  12,  13,  15,  18,  21  and  30  for  /?e  =  100,  d„ff  =  Q,  (T=1,  and 
v„ax  —  rjilna)  —  0.4.  The  contour  values  of  the  pseudo-streamlines  are  0,  ±0.02, 
±0.1  and  ±0.3.  The  contour  values  of  the  vorticity  are  ±0.4,  ±0.5,  ±0.8,  ±1.4  and 
±2,  with  the  highest  magnitude  at  the  sphere  surface.  The  solid  and  dotted  lines  in  the 
figures  represent  positive  and  negative  values.  Figures  6(a)-6(y)  show  the  pressure 
coefficient,  2{p—p^)/pUl„  and  shear  stress  coefficient,  Ir^JpUl^,  around  the  sphere  in 
the  principal  plane,  respectively,  at  /  =  1,  6,  9,  10,  11,  12,  13,  15,  18  and  21.  Note  that 
figure  6(a)  is  for  /  =  1  which  corresponds  to  figure  5(6). 

At  /  =  0,  figure  5  (a)  shows  that  a  spherical  vortex  sheet  is  generated  around  the 
sphere  owing  to  the  no-slip  condition  at  the  sphere  surface.  The  subsequent  figures 
show  that  the  vortex  sheet  around  the  sphere  is  advected  downstream  as  well  as 
diffused  outwards  from  the  sphere.  The  vorticity  on  the  edge  of  the  vortex  core  is  0.4 
at  f  =  0  for  =  2.51  and  a  =  1  which  correspond  to  =  rjiltia)  =  0.4.  The 
vortex  tube  is  initially  cylindrical  and  thus  should  appear  as  a  circle  in  the  principal 
plane.  However,  the  vortex  tube  in  figure  5(a)  is  not  an  exact  circle  because  the  grid 
resolution  is  relatively  coarse  at  the  initial  location  of  the  vortex  tube  which  is 
far  upstream  from  the  sphere  and  the  linear  interpolation  is  used  to  draw  the  contour 
lines.  However,  we  calculate  analytically  the  exact  velocity  and  pressure  fields  induced 
by  the  vortex  tube  by  using  equations  (1)  and  (2),  and  prescribe  them  as  initial 
conditions.  Therefore,  the  magnitudes  of  the  initial  velocity  components  at  a  given 
location  (x,  z)  are  fixed  no  matter  what  grid  distribution  is  used.  Thus,  the  circulation 
around  a  large  circle  enclosing  the  vortex  tube  remains  the  same  as  that  of  the  vortex 
tube.  The  velocity  and  pressure  fields  as  a  function  of  time  are  almost  not  affected  by 
the  initial  vortex  tube  shape  obtained  by  the  linear  interpolation.  The  line  connecting 
the  front  and  rear  stagnation  points  in  the  standard  axisymmetric  flow  over  a  single 
sphere,  which  is  the  .x  =  0  line  in  the  principal  plane,  will  be  used  as  a  reference  line. 
We  refer  to  the  region  above  the  line  as  ‘upper’  region  and  that  below  the  line  as 
‘lower’  region. 

For  0  <  r  <  9,  the  vortex  tube  is  upstream  of  the  sphere  as  shown  in  figures 
5  (6)-5  (</).  The  vortex  tube  rotating  counterclockwise  produces  downwash  upstream  of 
itself  and  upwash  downstream.  Therefore,  the  front  stagnation  point  on  the  sphere  is 
shifted  below  the  plane  x  =  0  owing  to  the  upwash,  and  thus,  the  fluid  particles  in  the 
upper  left-hand  region  move  faster  than  do  those  in  the  lower  left  region  of  the  sphere. 
As  a  consequence,  lower  pressure  and  higher  shear  stress  act  in  the  upper  left-hand 
region  compared  to  the  lower  left-hand  region  as  shown  in  figures  6(a)^(c),  and  this 
causes  a  positive  lift  force  on  the  sphere.  Note  that  in  figure  6  the  clockwise  direction 
is  considered  positive  for  the  shear  stress  in  the  upper  region  of  the  sphere,  and  the 
counterclockwise  direction  is  considered  positive  for  the  shear  stress  in  the  lower 
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region.  On  the  other  hand,  the  shift  of  the  front  stagnation  point  below  the  plane  a:  =  0 
causes  the  fluid  particles  to  continue  to  accelerate  after  6  =  90°  at  the  bottom  of  the 
sphere  but  to  begin  to  decelerate  before  6  =  90°  at  the  top  of  the  sphere  where  6  is 
measured  from  the  negative  z-axis,  as  shown  in  the  pressure  distribution  around  the 
sphere  in  the  principal  plane  in  figures  6(a)-6(c).  Thus  the  fluid  particles  move  faster 
in  the  bottom  and  lower  right-hand  regions  than  in-  the  top  and  upper  right-hand 
regions  of  the  sphere.  As  a  consequence,  lower  pressure  and  higher  shear  stress  act  in 
the  bottom  and  lower  right-hand  regions  compared  to  the  top  and  upper  right-hand 
regions  as  shown  in  figures  6(a)-6(c),  and  this  causes  the  fluid  particles  turning  around 
the  upper  eddy  to  be  pushed  into  the  lower  region  of  the  near  wake  as  shown  in  figures 
5(c)  and  5{d).  Figures  5(c)  and  5{d)  also  show  that  the  upper  eddy  is  formed  by  the 
fluid  separating  on  the  upper  portion  of  the  sphere  as  in  the  case  of  axisymmetric  flow 
past  a  sphere  without  the  presence  of  the  vortex  tube.  On  the  other  hand,  the  lower 
eddy  is  not  formed  by  the  fluid  separating  on  the  lower  portion  of  the  sphere,  but  rather 
by  the  fluid  turning  around  the  upper  eddy  and  being  entrained  by  the  lower  flow.  This 
lower  eddy  is  detached  from  the  sphere.  A  portion  of  the  fluid  moving  around  the  top 
of  the  sphere  passes  between  the  detached  lower  eddy  and  the  sphere.  A  similar  flow 
pattern  was  found  by  Kim,  Elghobashi  &  Sirignano  (1993)  in  their  study  of  three- 
dimensional  flow  over  two  spheres  placed  side  by  side. 

For  9  <  /  <  10.  the  figures  5(rf)  and  5(c)  shows  that  the  vortex  tube  contacts  the 
boundary  layer  of  the  sphere. 

For  10  <  t  <  13,  the  figures  5(/)-5(/»)  show  that  the  vortex  tube  goes  around  the 
bottom  of  the  sphere.  The  vortex  tube  is  now  downstream  of  the  front  stagnation  point 
in  the  axisymmetric  flow  past  a  sphere  and  produces  downwash  on  the  sphere. 
Therefore,  the  front  stagnation  point  on  the  sphere  is  shifted  above  the  plane  x  =  0, 
and  thus  the  fluid  particles  in  the  lower  left-hand  region  move  faster  than  do  those  in 
the  upper  left-hand  region  of  the  sphere.  As  a  consequence,  lower  pressure  and  higher 
shear  stress  act  in  the  lower  left-hand  region  compared  to  the  upper  left-hand  region 
as  shown  in  flgures  6(c)-6(g),  and  this  causes  the  negative  lift  force  on  the  sphere.  On 
the  other  hand,  the  shift  of  the  front  stagnation  point  above  the  plane  x  =  0  causes  the 
fluid  particles  to  move  faster  in  the  boundary  layer  of  the  top  and  upper  right-hand 
regions  compared  to  that  of  the  bottom  and  lower  right-hand  regions  of  the  sphere  as 
shown  in  figure  7,  which  shows  the  tangential  velocity  profiles,  u„{r),  at  ^  =  90°  on  the 
top  and  the  bottom  of  the  sphere  in  the  principal  plane  at  r  =  12,  and  this  causes  the 
higher  shear  stress  in  the  top  and  upper  right-hand  regions  compared  to  the  bottom 
and  lower  right-hand  regions  as  shown  in  the  figures  6  if)  and  6  (g).  However,  during 
this  time  period,  the  pressure  distributions  at  the  top  and  bottom  of  the  sphere  in  the 
principal  plane  have  different  features  from  those  of  the  shear  stress  owing  to  the 
following  reason.  The  counterclockwise  vortex  tube  in  the  uniform  stream  produces  a 
flow  field  in  which  the  fluid  velocity  is  less  than  that  of  the  uniform  stream  above  the 
vortex  tube  and  higher  than  that  of  the  uniform  stream  below  the  vortex  tube  with 
respect  to  fixed  coordinate  system  in  space.  Owing  to  this  shear  flow,  the  fluid  velocity 
on  the  edge  of  (and  outside)  the  boundary  layer  at  the  bottom  of  the  sphere  is  larger 
than  that  at  the  top  of  the  sphere  as  shown  in  figure  7,  and  thus  the  pressure  at  the 
bottom  of  the  sphere  is  lower  than  that  at  the  top  of  the  sphere.  This  pressure  difference 
causes  the  fluid  particles  turning  around  the  top  of  the  sphere  to  be  pushed  into  the 
lower  region  of  the  wake  forming  an  S-shaped  path  (figures  5(g)  and  5(A)).  The 
combined  effect  of  the  upward  shift  of  the  front  stagnation  point  owing  to  downwash 
of  the  vortex  tube  and  the  velocity  difference  between  the  top  and  the  bottom  of  the 
sphere  owing  to  the  shear  flow  induced  by  the  vortex  tube  results  in  a  higher  magnitude 
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Figure  7.  Tangential  velocity  profile,  ujir),  at  ^  =  90®  on  the  top  and  the  bottom 
of  the  sphere  in  the  principal  plane  at  /  =  12. 


of  the  maximum  negative  force  than  that  of  the  maximum  positive  force,  as  will  be 
shown  in  detail  in  §3.2.2.  The  upper  separation  eddy  becomes  smaller  during  this  time 
period,  because  the  pressure  difference  between  the  upper  and  lower  wake  just 
downstream  of  the  sphere  becomes  larger  when  the  vortex  tube  passes  the  plane  z  =  1 
(tangent  to  the  rear  stagnation  point  in  axisymmetric  flow),  and  more  fluid  particles  are 
pushed  into  the  lower  wake  just  behind  the  sphere.  At  t  =  12  and  13,  no  separation 
eddies  appear  in  the  wake  as  shown  in  figures  5(g)  and  5(A),  and  the  flow  does  not 
separate  in  the  upper  region  of  the  sphere  as  shown  in  figures  5(g),  5(A),  6(/)  and  6(g). 
We  note  that  the  separation  point  in  the  principal  plane  is  the  point  at  which  the  shear 
stress  vanishes. 

The  reason  for  the  passage  of  the  counterclockwise-rotating  vortex  tube  around  the 
bottom  of  the  sphere  rather  than  around  the  top  is  as  follows.  First,  note  that  the  well- 
known  two-dimensional,  inviscid  case  of  a  vortex  interacting  with  a  cylinder  has  a 
counterclockwise  (clockwise)  rotating  vortex  advecting  clockwise  (counterclockwise) 
around  the  cylinder.  In  our  case,  the  opposite  behaviour  suggests  that  viscosity  is 
important  in  this  phenomenon.  Note  further  that  the  vorticity  levels  associated  with 
the  viscous  boundary  layer  on  the  sphere  are  greater  than  those  associated  with  the 
tube.  When  the  counterclockwise-rotating  vortex  tube  comes  close  to  the  sphere 
boundary  layer,  it  augments  the  magnitude  of  the  edge  velocity  in  the  lower  boundary 
layer  and  reduces  the  edge  velocity  in  the  upper  boundary  layer.  The  result  is  a  higher 
strength  vorticity  in  the  lower  boundary  layer  than  in  the  upper  boundary  layer  (see  the 
vorticity  contours  in  figures  5(d)-5  (/)).  (The  magnitude  of  the  highest  vorticity  in  the 
lower  boundary  layer  is  1 5  %  higher  than  that  in  the  upper  boundary  layer  at  t  =  9.) 
Consequently,  the  vorticity  in  the  lower  boundary  layer  induces  a  velocity  in  the 
downward  direction  at  the  location  of  the  vortex  tube  with  higher  magnitude  than  that 
induced  by  the  vorticity  in  the  upper  boundary  layer.  This  downward  induced  velocity 
advects  the  vortex  tube  below  the  sphere. 

For  13  <  t  <  19,  the  vortex  tube  is  downstream  of  the  sphere  as  shown  in  figures  5(0 
and  5(j)  and  produces  downwash  on  the  sphere.  Therefore,  the  negative  lift  force  acts 
on  the  sphere  owing  to  the  shift  of  the  front  stagnation  point  above  the  plane  x  =  0  in 
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(b) 


Figure  8.  Two  views  of  a  three-dimensional  contour  surface  of  =  0.2  at  t  =  20  for  the  flow 
depicted  in  figure  5 ;  (n)  a  side  view  looking  normal  to  the  principal  plane  (b)  a  view  looking  down 
with  an  acute  angle  toward  the  (3>,z)-plane. 

a  similar  manner  as  for  10  <  /  <  13,  but  the  negative  lift  force  is  reduced  as  the  vortex 
tube  moves  further  downstream.  The  shift  of  the  front  stagnation  point  above  the 
plane  x  =  0  causes  the  fluid  particles  to  move  faster  in  the  top  and  upper  right-hand 
regions  than  in  the  bottom  and  lower  right-hand  regions  of  the  sphere,  and  this  causes 
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the  lower  pressure  and  higher  shear  stress  in  the  top  and  upper  right-hand  regions 
compared  to  those  in  the  bottom  and  lower  right-hand  regions  as  shown  in  figures  6{h) 
and  6(/).  However,  because  the  vortex  tube  is  still  intersecting  the  near  wake  and  thus 
producing  strong  downwash  in  the  near  wake,  the  fluid  particles  turning  around  the 
top  of  the  sphere  are  pushed  into  the  lower  region  of  the  near  wake.  This  allows  no 
room  for  the  lower  eddy  to  grow.  On  the  other  hand,  the  upper  eddy  grows  as  the 
vortex  tube  moves  downstream  because  the  fluid  particles  turning  around  the  top  of  the 
sphere  experience  less  force  pushing  them  into  the  lower  region  of  the  near  wake. 

For  t  ^  20,  the  vortex  tube  is  far  downstream  of  the  sphere  as  shown  in  figures  5(k) 
and  5(/)  and  produces  weak  downwash  on  the  sphere,  and  thus  the  lift  force  on  the 
sphere  is  almost  zero  as  will  be  shown  in  §3.2.2.  The  weak  downwash  causes  the  front 
stagnation  point  on  the  sphere  to  be  shifted  slightly  above  the  plane  .v  =  0,  and  thus 
the  fluid  particles  move  slightly  faster  in  the  top  and  upper  right-hand  regions  than  in 
the  bottom  and  lower  right-hand  regions  of  the  sphere,  and  this  causes  the  higher  shear 
stress  and  lower  pressure  in  the  top  and  upper  right-hand  regions  as  shown  in  figure 
6(j).  Now,  the  downward  force  due  to  the  vortex  tube  is  very  weak  in  the  near  wake 
because  the  vortex  tube  is  far  downstream.  Therefore,  the  lower  eddy  grows,  and  owing 
to  the  lower  pressure  in  the  upper  region  of  the  sphere,  the  fluid  particles  turning 
around  the  lower  eddy  are  pushed  up  into  the  upper  near  wake  as  shown  in  figures  5ik) 
and  5(1). 

We  now  examine  a  three-dimensional  view  of  the  vortex  tube  by  considering  the  y- 
component  of  vorticity  vector.  Figures  8(a)  and  i{b)  show  two  views  of  a  three- 
dimensional  contour  surface  of  Wj,  =  0.2  at  /  =  20  for  the  flow  depicted  in  figure  5. 
Figure  8(a)  shows  a  side  view  looking  normal  to  the  principal  plane,  whereas  figure 
8  (b)  shows  a  view  looking  down  with  an  acute  angle  toward  the  (y,  z)-plane.  The  ellipse 
in  figure  Sib)  is  the  boundary  of  the  spherical  computational  domain  viewed  at  an  angle. 
It  appears  as  a  circle  when  viewed  normal  to  the  principal  plane.  The  sphere  is  at  the 
centre  of  the  domain  in  figures  8(a)  and  S(b).  Figure  8(6)  shows  that  the  portion  of  the 
vortex  tube  in  the  principal  plane  is  retarded,  owing  to  its  interaction  with  the  sphere, 
compared  with  the  rest  of  the  vortex  tube  (in  the  (>',z)-plane  with  its  axis  parallel  to 
the  y-axis)  outside  the  principal  plane.  By  measuring  the  radial  extent  of  the  contour 
surface  (of  (Oy  =  0.2  at  t  =  20)  in  figures  8(a)  and  8(6),  we  find  that  the  maximum 
radius  of  the  contour  surface  outside  the  principal  plane  is  1.58  which  is  very  close  to 
the  value  of  1.6  taken  from  figure  3(6). 

3.2.2.  Lift,  moment,  and  drag  coefficients  and  effect  of  tube  circulation 

Figure  9  shows  the  lift  coefficients  of  the  sphere  as  a  function  of  time  for  Re  =  100, 
=  0,  and  o-  =  1 .  The  lift  coefficients  are  computed  for  four  different  maximum 
fluctuation  velocities  iv^„^  =  r,,/(27tcr))  due  to  the  vortex  tube,  with  magnitudes  equal 
to  0.1,  0.2,  0.3  and  0.4  (normalized  by  free-stream  velocity).  Owing  to  the  sudden 
placement  of  the  sphere  into  the  stream,  it  initially  takes  a  small  time  (0  <  /  <  0.6)  for 
the  initial  flow  perturbations  to  vanish. 

As  discussed  earlier,  when  the  vortex  tube  approaches  the  sphere  (0  <  /  ^  9.4),  it 
produces  upwash  resulting  in  a  positive  lift  force  on  the  sphere.  The  maximum  positive 
lift  coefficient  Q  occurs  at  t  =  7.2.  On  the  other  hand,  when  the  vortex  tube 
passes  the  sphere,  it  produces  downwash  and  high  fluid  velocity  near  the  bottom  of  the 
sphere  resulting  in  a  negative  lift  force.  The  magnitude  of  the  negative  lift  is  greater 
than  the  positive  lift  (figure  9).  The  maximum  negative  lift  coefficient  Q occurs 
at  /=  1 1.8  about  when  the  centre  of  the  vortex  tube  passes  the  plane  z  =  1.  The  lift 
coefficient  is  linearly  proportional  to  the  maximum  fluctuation  velocity  (or  the 
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Figure  9.  Lift  coefficients  of  the  sphere  as  a  function  of  time  and  for 
Re  =  100,  =  0  and  (r  =  1 . 


circulation  of  the  vortex  tube  for  constant  vortex  core  radius)  until  the  vortex  tube 
contacts  the  sphere  boundary  layer  (t  ^  9.4).  The  maximum  positive  lift  coefficient 
Q. maxi  is  expressed  by 

^L.maxl  ~  ^max^  (18) 

where  the  proportionality  constant  c  =  0.8.  For  t  >  9.4,  the  relation  between  the  lift 
coefficient  and  v^^x  deviates  slightly  from  linearity,  but  the  maximum  negative 
lift  coefficient  Q  ^0^2  is  linearly  proportional  to  v^^x  "'•th  c  =  —  1.66.  After  the  lift 
coefficient  reaches  its  maximum  negative  value,  it  decays  quickly  towards  zero  because 
the  vortex  tube  vorticity  is  diffused  in  the  sphere  wake.  The  time  averaged  lift 
coefficient  (averaged  over  a  time  span  between  t  =  0.6  and  the  maximum  time  24)  for 
all  values  of  is  small  and  negative  ((9(10"*)).  As  mentioned  earlier,  the  behaviour 
of  Q(0  during  the  period  0  <  r  <  0.6  is  influenced  by  the  initial  flow  perturbation,  and 
thus  its  value  during  this  initial  period  is  excluded  from  the  averaging  process.  The  root 
mean  square  Q.rm,  of  the  lift  coefficient  as  a  function  of  time  is  also  linearly 
proportional  to  v^ax  with  c  =  0.65  as  will  be  shown  in  table  4. 

Figure  10  shows  the  temporal  development  of  the  moment  coefficients  for  the  sphere 
under  the  same  conditions  as  figure  9.  The  moment  coefficients  are  obtained  for  four 
different  values  of  =  0.1,  0.2,  0.3  and  0.4. 

As  the  vortex  tube  approaches  the  sphere,  the  downward  shift  of  the  front  stagnation 
point  (due  to  the  upwash)  causes  higher  shear  stress  in  the  upper  left-hand  region 
compared  to  the  lower  left-hand  region  generating  a  negative  (clockwise)  torque.  At 
the  same  time,  the  downward  shift  causes  higher  shear  stress  in  the  bottom  and  lower 
right-hand  regions  compared  to  the  top  and  upper  right-hand  regions  as  explained  in 
§3.2.1  generating  a  positive  torque.  The  two  torques  compete  with  each  other  and 
result  in  a  net  weak  torque  in  the  interval  0  <  (  <  9. 

As  the  vortex  tube  passes  the  sphere  (9  <  t  ^  14),  the  upward  shift  of  the  front 
stagnation  point  (due  to  the  downwash)  causes  higher  shear  stress  in  the  lower  left- 
hand  region  compared  to  the  upper  left-hand  region  generating  a  positive  torque.  At 
the  same  time,  the  upward  shift  caused  higher  shear  stress  in  the  top  and  upper  right- 
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Figure  10.  Moment  coefficients  of  the  sphere  under  the  same  conditions  as  figure  9. 


hand  regions  compared  to  the  bottom  and  lower  right-hand  regions  as  explained  in 
§3.2.1  generating  a  negative  torque.  However,  the  effect  of  this  negative  torque  is 
diminished  by  the  shear  flow  induced  by  the  vortex  across  the  sphere  which  produces 
high  shear  stress  at  the  bottom  of  the  sphere.  As  a  consequence,  a  net  high  positive 
torque  acts  on  the  sphere.  The  maximum  positive  moment  coefficient  C„  occurs 
at  /=  11.4.  C„  is  linearly  proportional  to  with  a  proportionality  constant 
c  =  0.14. 

When  the  vortex  tube  is  relatively  far  downstream  from  the  sphere  (t  >  15),  the 
positive  torque  due  to  the  shear  stress  in  the  lower  left-hand  region  competes  with  the 
negative  torque  due  to  the  shear  stress  in  the  top  and  upper  right-hand  regions.  This 
results  in  a  net  weak  negative  torque  which  becomes  smaller  as  the  vortex  tube  moves 
farther  downstream.  We  note  that  the  torque  depends  only  on  the  distribution  of  the 
shear  stresses  and  is  relatively  small  compared  to  the  lift  force. 

The  time-averaged  moment  coefficient  (averaged  over  a  time  span  between  t  =  0.6 
and  24)  for  all  values  of  is  nearly  zero  (0(10”*)),  and  the  root  mean  square  C,,, 
of  the  moment  coefficient  is  approximately  linearly  proportional  to  with  c  =  0.05. 

Figure  1 1  shows  the  drag  coefficients  of  the  sphere  as  a  function  of  time  for  the  same 
conditions  as  figure  9.  The  drag  coefficients  are  computed  for  four  different  values  of 
t;„„^  =  0.1,  0.2,  0.3  and  0.4. 

As  discussed  earlier,  the  sudden  placement  of  the  sphere  in  the  flow  results  in  initially 
large  values  of  shear  stress  and  pressure  on  the  sphere,  and  hence  a  large  drag  as  shown 
in  figure  11.  Figure  5(e)  shows  that  at  about  t=  10  the  centre  of  the  vortex  tube  is 
located  near  the  front  stagnation  point  which  is  slightly  below  the  point  (.Y,y,z)  = 
(0,0,  - 1).  Owing  to  the  low  pressure  at  the  centre  of  the  vortex  tube,  the  pressure 
coefficient  at  the  front  stagnation  point  (Cp„  =  0.818)  is  lower  than  that  of  the 
axisymmetric  flow  past  a  sphere  without  the  vortex  tube  (Cp„  =  1.107)  as  shown  in 
figure  6 (if).  Also,  the  maximum  shear  stresses  in  the  upper  and  lower  regions  are  lower 
than  that  of  the  axisymmetric  flow  without  the  vortex  tube.  This  causes  the  drag  on  the 
sphere  to  be  lower  than  that  of  the  axisymmetric  flow  without  the  vortex  tube.  As  the 
vortex  tube  moves  around  the  bottom  of  the  sphere,  the  front  stagnation  point  is 
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Figure  11.  Drag  coefficients  of  the  sphere  under  the  same  conditions  as  figure  9. 

shifted  above  the  plane  .v  =  0,  owing  to  the  downwash.  Consequently,  high  pressure 
and  higher  shear  stress  act  in  the  upper  and  lower  left-hand  regions,  respectively,  as 
explained  earlier  in  §3. 2.1.  This  increases  the  drag  during  the  periods  10  <  r  <  13.4.  For 
t>  13.4,  the  drag  approaches  that  of  the  axisymmetric  Row  as  the  vortex  tube  moves 
further  downstream. 

The  time-averaged  value  of  the  deviation  of  the  drag  coefficient  from  that  of  the 
axisymmetric  flow  past  a  sphere  for  all  values  of  is  nearly  zero  (0(1 0”“')).  The 
unsteady  drag  coefficient  of  the  axisymmetric  flow  past  a  sphere  was  computed  for  a 
sphere  suddenly  placed  in  the  uniform  stream  without  the  vortex  tube. 

3.3.  Effects  of  the  offset  distance 

We  examine  the  effects  of  the  offset  distance  on  the  flow  field  by  varying  d„ff  while  using 
the  same  flow  conditions  as  in  §3.2. 

3.3.1.  Offset  distance  1  ^  d„ff  ^4 

The  temporal  behaviours  of  the  lift  and  moment  coefficients  of  the  sphere  for 
d„ff  =  1  are  similar  to  those  in  the  case  of  d„ff  =  0.  The  main  feature  distinguishing  the 
case  of  d^ff  =  1  from  that  of  d„ff  =  0  is  that  in  the  former,  the  vortex  tube  splits  into 
two  parts  when  the  vortex  tube  passes  the  sphere.  The  attraction  of  the  vortex  tube  to 
the  positive  vorticity  in  the  boundary  layer  at  the  bottom  of  the  sphere  causes  some 
portion  of  the  vortex  tube  to  move  around  the  bottom  of  the  sphere,  whereas  the  other 
portion  moves  on  the  top  of  the  sphere,  as  shown  in  figure  12.  Figures  12(a)-12(/i) 
display  the  contour  lines  of  ,v-component  vorticity  in  the  principal  plane  at  /  =  9,  10, 

11,  12,  13,  15,  18  and  21  for  Re  =  100,  doff=  \.  a=  I,  and  =  0.4.  The  contour 

values  of  the  vorticity  are  ±0.4,  ±0.5,  ±0.8,  ±1.4  and  ±2.  Owing  to  its  longer 

interaction  with  the  sphere  for  =  I  than  for  d„ff  =  0,  the  magnitudes  of  the  lift  and 

moment  coefficients  of  the  sphere  are  close  to  those  in  the  case  of  d„ff  =  0  despite  its 

positive  offset  distance  initially.  Equation  (18)  is  approximately  valid  for  Q  „,oxp 

Cl.  max 2-  Q.  rm*-  Q/,  max^  and  C,,.  with  the  same  proportionality  constants  as  in  the  ‘  ' 

case  of  d„ff  =  0. 
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Figure  12.  Contour  lines  of  ji'-component  vorticity  in  the  principal  plane  at  (a)  t  =  9,  (b)  10, 

(c)  1 1,  (d)  12,  (e)  13,  (/)  15,  (g)  18  and  (h)  21  for  Re  =  100,  =  1,  a  =  1  and  =  0.4. 

Figure  13(a)  shows  the  drag  coefficients  of  the  sphere  as  a  function  of  time  for 
Re  =  100.  dgff  =  1,  and  cr  =  1.  The  drag  coefficients  are  obtained  with  two  different 
maximum  fluctuation  velocities  due  to  the  vortex  tube,  £)„,„j  =  0.1  and  0.2.  The 
temporal  behaviour  of  the  drag  coefficients  is  different  from  that  of  the  case  of  d„ff  =  0. 
The  time-averaged  value  of  the  deviation  of  the  drag  coefficient  from  that  of  the 
axisymmetric  flow  past  a  sphere  for  all  values  of  nearly  zero  but  increases 

linearly  as  v^ax  increases.  The  time-averaged  drag  coefficient  ^ve  is  expressed  by 

^D.ave  ^ D,axi~^  (19) 

where  the  constant  p  =  0.2,  and  is  the  time-averaged  value  of  the  drag  coefficient 
in  the  case  of  axisymmetric  flow  {v^^x  =  0)-  The  drag  coefficients  reach  their  maximum 
at  about  t  =  10.  The  maximum  drag  coefficient  i®  expressed  also  by  (19)  but 

with  /?  =  0.62,  and  ^xt  here  is  the  local  value  of  the  axisymmetric  drag  coefficient  at 
the  time  of  At  about  r  =  10,  the  centre  of  the  vortex  tube  is  located  above  the 

front  stagnation  point.  Thus,  the  induced  velocity  due  to  the  vortex  tube  adds  its 
magnitude  to  the  base  flow  along  the  stagnation  streamline,  and  so  the  dynamic 
pressure  ahead  of  the  front  stagnation  point  becomes  higher  than  that  of  the 
axisymmetric  flow  past  a  sphere.  This  causes  the  pressure  at  the  stagnation  point  and 
the  shear  stresses  in  the  upper  and  lower  left-hand  regions  to  be  higher  than  those  of 
the  axisymmetric  flow  past  a  sphere.  As  a  consequence,  the  drag  is  increased.  When  the 
offset  distance  is  negative,  the  reverse  phenomena  would  occur,  and  the  drag  would  be 
decreased.  This  will  be  discussed  in  §3.3.2. 

Figures  14(a)  and  14(6)  display  the  contour  lines  of  y-component  vorticity  in  the 
principal  plane  at  /  =  9  and  12  for  Re  =  100,  d„ff  =  2,  a  =  1,  and  v^^x  =  O  The 
contour  values  of  vorticities  are  the  same  as  those  of  previous  sections.  Figures  14(a)  and 


Figure  13.  Drag  coefficients  of  the  sphere  as  a  function  of  time  and  for  {a) 

(c)  3  and  (d)  4  with  Re  =  100  and  <r  =  1. - ,  axisymmetric; - ,  =  0.1 ; - , 


Figure  14.  Contour  lines  of  y-component  vorticity  in  the  principal  plane  at 
(a)  t  =  9  and  (b)  12  for  Re  =  100,  d.„  =  2.  tr  =  1  and  =  0.4. 
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Figure  15.  Magnitude  of  and  C„  as  a  function  of  for 

Re  =  100,  =  1  and  =  0.2. 


14(6)  show  that  the  vortex  tube  passes  above  the  sphere.  The  behaviour  of  the  lift 
coefficients  with  time  for  d^ff  =  2,  3  and  4  is  similar  to  that  of  the  case  of  d^ff  =  0  and 
1 .  However,  their  magnitudes  are  smaller  than  those  for  dgff  =  0  and  1  and  decay  with 
d^ff  exponentially  as  shown  in  figure  15,  where  the  magnitude  of  the  negative  maximum 
lift  coefficient  and  the  maximum  moment  coefficient  for  =  100,  cr=l,  and 
^max  =  0-2  are  presented  as  a  function  of  d^ff.  The  positive  maximum  lift  coefficient 
Cl. maxi  for  doff  >2  is  expressed  by 

Cl.  max  1  =  <^1  fmax  eXP  (C^  WoffW  (20) 

where  q  =  0.99  and  Cj  = -0.3.  The  negative  maximum  lift  coefficient  Cl. maxi  ‘S 
expressed  by  (20)  with  =  -2.64  and  =  -0.38,  and  the  r.m.s.  lift  coefficient 
is  expressed  also  by  (20)  with  Cj  =  0.88  and  =  —0.28. 

The  behaviour  of  the  moment  coefficients  with  time  for  d„ff  =  2, 3  and  4  is  also 
similar  to  that  of  the  case  of  doff  =  0  and  1 .  However,  their  magnitudes  are  smaller  than 
those  for  doff  =  0  and  1  and  decay  with  a  negative  power  of  dgff  as  shown  in  figure  15. 
The  maximum  moment  coefficient  Cnj  ,„ox  for  dgff  ^  2  is  expressed  by 

Cm. max  ^i^max\doff\  ’  (21) 

where  =  0.185  and  m=  — 1.501.  The  r.m.s.  moment  coefficient  ^ms  is  expressed 
by  (21)  with  Cg  =  0.056  and  m=  —  1.185. 

Figures  13(6),  13(c)  and  13(</)  show  the  drag  coefficients  of  the  sphere  as  a  function 
of  time  for  doff  =  2,  3  and  4,  respectively,  with  Re  =  100  and  cr=  1.  The  drag 
coefficients  are  obtained  with  two  different  maximum  fluctuation  velocities  due  to  the 
vortex  tube,  v^^x  =  01  ^nd  0.2.  The  drag  coefficients  reach  their  maximum  at  about 
t  =  10.  The  maximum  drag  coefficient  C^  for  d„ff  =  2  is  higher  than  that  for 
doff  =  1  because  the  magnitude  of  the  induced  velocity  added  to  the  base  flow  along  the 
stagnation  streamline  for  =  2  is  higher  than  that  for  dgff  =  1 .  We  note  that  the 
radius  of  the  vortex  core  is  greater  than  unity  at  t  =  10  owing  to  the  diffusion  (and  the 
maximum  induced  velocity  occurs  at  the  edge  of  the  vortex  core)  as  shown  in  figure  3  (a) 
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for  Re  =  100,  <r  =  1,  and  =  0.4.  Cj,  decreases  as  d„ff  increases  above  2.  C„ 
for  d„ff  =  2,3  and  4  are  expressed  by  (19)  with  /?  =  0.9,  0.8  and  0.68,  respectively. 
Cy,  ave  for  ^off  =  2, 3  and  4  are  expressed  also  by  (19)  with  /^  =  0.31,  0.33  and  0.32, 
respectively.  We  note  that  the  magnitude  of  the  deviation  of  the  drag  coefficient  from 
that  of  the  axisymmetric  flow  decays  slowly  with  d^fp  in  contrast  with  fast  decay  of  the 
lift  and  moment  coefficients  with  d^jf. 

3.3.2.  Offset  distance  —  1  ^  d^ff  ^  —4 

Note  that  the  sign  reversal  of  the  initial  tube  vorticity  with  the  offset  distance  kept 
positive  is  a  mirror  image  of  the  case  where  the  sign  of  the  offset  distance  is  changed 
and  the  sign  of  the  initial  vorticity  is  kept  constant.  Therefore,  we  consider  only  change 
in  sign  of  the  offset  distance  and  keep  the  counterclockwise  rotation. 

The  behaviour  of  the  lift  coefficients  with  time  for  —  1  ^  dgff  ^  —  4  is  similar  to  that 
of  the  case  of  d„ff  =  0.  However,  their  magnitudes  are  smaller  than  that  for  d^ff  =  0 
with  the  same  ^nd  decay  exponentially  with  d^ff  as  shown  in  figure  15.  The 
positive  maximum  lift  coefficient,  the  negative  maximum  lift  coefficient,  and  the  r.m.s. 
lift  coefficient  for  d„ff  <  - 1  are  expressed  by  (20)  with  c,  =  0.942  and  —  —0.295, 
Cl  =  -  1.95  and  Cj  =  -0.35,  and  Cj  =  0.74  and  Cj  =  -0.27,  respectively. 

The  behaviour  of  the  moment  coefficients  with  time  for  —  1  ^  d„ff  ^  —  4  is  similar  to 
that  of  the  case  of  d^j^  =  0.  However,  their  magnitudes  are  smaller  than  that  for 
d^ff  =  0  with  the  same  and  decay  with  a  negative  power  of  \d„ff\  as  shown  in 
figure  15.  The  maximum  moment  coefficient  and  the  r.m.s.  moment  coefficient  are 
expressed  by  (21)  with  Cj  =  0.09  and  m  =  — 1.264,  and  =  0.0318  and  m  =  — 1.047, 
respectively. 

Figure  16(a)  shows  the  drag  coefficients  of  the  sphere  as  a  function  of  time  for 
Re  =  100,  d„ff  =  - 1,  and  cr  =  1.  The  drag  coefficients  are  obtained  with  two  different 
maximum  fluctuation  velocities  due  to  the  vortex  tube,  c„„^  =  0.1  and  0.2.  The 
behaviour  of  the  drag  coefficients  with  time  is  different  from  that  of  the  case  of  d„ff  =  0. 
The  time-averaged  value  of  the  deviation  of  the  drag  coefficient  from  that  of  the 
axisymmetric  flow  past  a  sphere  at  each  maximum  fluctuation  velocity  is  not  near  zero 
but  is  increased  linearly  in  negative  value  as  the  maximum  fluctuation  velocity  becomes 
higher.  The  minimum  drag  coefficients  occur  at  about  /  =  10.  The  minimum  drag 
coefficient  and  the  time-averaged  drag  coefficient  are  expressed  by  (19)  with  p  =  -0.78 
and  -0.2,  respectively. 

At  about  /  =  10,  the  centre  of  the  vortex  tube  is  located  below  the  front  stagnation 
point.  Thus,  the  induced  velocity  due  to  the  vortex  tube  subtracts  its  magnitude  from 
the  base  flow  along  the  stagnation  streamline,  and  so  the  dynamic  pressure  ahead  of 
the  front  stagnation  point  becomes  lower  than  that  of  the  axisymmetric  flow  past  a 
sphere.  This  causes  the  pressure  at  the  front  stagnation  point  and  the  shear  stresses  in 
the  upper  and  lower  left  regions  to  be  lower  than  those  of  the  axisymmetric  flow  past 
a  sphere.  As  a  consequence,  the  drag  is  decreased.  From  this  result,  we  deduce  that  if 
the  sphere  were  free  to  move  rather  than  fixed,  it  would  experience  lower  drag  than 
that  of  a  sphere  subjected  to  an  axisymmetric  flow  unless  the  initial  offset  distance  is 
large  positive.  The  lower  drag  will  be  caused  by  the  upward  motion  of  the  sphere  owing 
to  the  upwash  when  the  vortex  tube  approaches  it,  and  thus  the  centre  of  the  vortex 
tube  will  be  located  below  the  front  stagnation  point  of  the  sphere.  This  will  cause 
lower  dynamic  pressure  ahead  of  the  front  stagnation  point. 

Figures  16(6)-16(i/)  show  the  drag  coefficients  of  the  sphere  as  a  function  of  time  for 
dgff  =  —  2,  —  3  and  —  4  with  Re  =  1 00  and  er  =  1 .  The  minimum  drag  coefficients  occur 
at  about  /  =  10.  The  minimum  drag  coefficient  Q,  for  d^ff  =  -  2  is  lower  than  that 
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(a)  (c) 


Figure  16.  Drag  coefficients  of  the  sphere  as  a  function  of  time  and  for  (a)  =  - 1.  (b) 

-2,  (c)  -3  and  (d)  -4  with  =  100  and  <r=l.  - .  axisymmetric; - .  =  0.1; 

■  •  •  •.  =  0.2. 


for  d„ff  =  - 1,  because  the  magnitude  of  the  induced  velocity  subtracted  from  the  base 
flow  along  the  stagnation  streamline  is  higher  for  =  -  2  compared  to  d„ff  =  —  1. 
The  magnitude  of  increases  as  d„jf  decreases  below  —2.  The  minimum  drag 

coefficients  for  d^jj  =  —2,  —3  and  —4  are  expressed  by  (19)  with  /?  =  —0.9,  —0.8  and 
—0.68,  respectively.  The  time-averaged  drag  coefficients  for  d^ff  =  -2,  -3  and  —4  are 
expressed  also  by  (19)  with  /?  =  —0.28,  —0.3  and  -0.29,  respectively.  We  note  that  the 
magnitude  of  the  deviation  of  the  drag  coefficient  from  that  of  the  axisymmetric  flow 
decays  slowly  with  d^ff,  in  contrast  with  fast  decay  of  the  lift  and  moment  coefficients 
with  d„ff. 

3.4.  Effects  of  the  size  of  the  vortex  tube 

We  examine  the  effects  of  the  size  of  the  vortex  tube  on  the  flow  field  by  performing 
computations  similar  to  those  in  §3.2  for  Re  =  100,  d„ff  =  0,  and  five  different  sizes  of 
the  vortex  tube,  0.25  <  cr  <  4  in  addition  to  the  base  case  cr  =  1 .  Each  simulation  is 
performed  with  two  different  values  of  =  0.1  and  0.3. 

Table  4  shows  the  maximum  positive  lift  coefficient,  the  maximum  negative  lift 
coefficient,  the  r.m.s.  lift  coefficient,  the  maximum  moment  coefficient,  and  the  r.m.s. 
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^  max 

c 

max  1 

^  I.,  max 2 

c 

rma 

(7  =  4 

^ .\I .  max 

Cu.xm- 

0.3 

0.337 

-0.566 

0.307 

0.036  1 

0.0161 

0.1 

0.113 

-0.190 

0.102 

(T  —  3 

0.0121 

0.00546 

0.3 

0.338 

-0.604 

0.299 

0.0414 

0.0169 

O.l 

0.114 

-0.201 

0.0997 

(T  =  2 

0.0139 

0.00577 

0.3 

0.318 

-0.611 

0.271 

0.0455 

0.0169 

0,1 

0.107 

-0.203 

0.0907 

0.0155 

0.00593 

0.3 

0.238 

-0.500 

(T  =  1 

0.195 

0.0410 

0.0138 

0.1 

0.0810 

-0.165 

0.0658 

<7  =  0.5 

0,0145 

0.00514 

0.3 

0.147 

-0.322 

0.119 

0.0290 

0.00953 

0.1 

0.0499 

-0.106 

0.0402 

tr  =  0.25 

0.0103 

0.00355 

0.3 

0.0815 

-0.177 

0.0660 

0.0175 

0.00587 

0.1 

0.0275 

-0.0583 

0.022 1 

0.00604 

0.00210 

Table  4.  Maximum  positive  lift  coefficient,  maximum  negative  lift  coefficient,  root  mean  square  of 
the  lift  coefficient,  maximum  moment  coefficient,  and  root  mean  square  of  the  moment  coefficient 
as  a  function  of  for  six  different  radii  of  the  vortex  tube,  (x  =  4,  3,  2,  1,  0.5  and  0.25. 


moment  coefficient  as  a  function  of  for  six  different  initial  radii  of  the  vortex  tube, 
(7  =  4,  3,  2,  1, 0,5  and  0.25.  All  the  coefficients  are  linearly  proportional  to  at  each 
(T.  When  cr  ^  2,  Q  and  C,^  become  independent  of  a,  but  the  magnitudes  of 
Q.  m<ix2’  C,/.  max'  Q/.  rm»  for  O'  =  4  are  Smaller  than  those  for  o'  =  2  and  3.  When 
O'  approaches  zero,  all  the  coefficients  tend  to  be  proportional  to  which  is  the 

circulation  of  the  vortex  tube  divided  by  2jt  =  27to'o„„_j).  For  example,  Q  ,.m« 
expressed  by 

Q.rm«  =  C,J^„ax  (2  ^  O' <  4), 

=  C2«^a,axO-"  (0.25  <  O' <2,  0.75^0^0.5),  (22) 

where  the  constant  c,  =  1  and  c._,  =  0.65,  and  n  depends  on  ir  and  should  approach 
unity  as  a  reaches  zero.  For  q  =  1.1  and  c,  =  0.8.  Q  and 

rm»  for  O'  ^  3  are  also  expressed  by  (22)  with  c,  =  —  2  and  =  — 1.66,  c,  =  0.14 
and  c.,  =  0.14,  and  c,  =  0.055  and  c.^  =  0.05,  respectively.  The  time-averaged  value  of 
the  deviation  of  the  drag  coefficient  from  that  of  the  axisymmetric  flow  past  a  sphere 
for  all  values  of  a  is  nearly  zero  (0(10'''')). 

We  note  that  Cy,  ,nax^  and  for  o-  =  4  are,  respectively,  smaller  than 

those  for  O'  =  2  and  3,  and  the  reason  is  explained  as  follows.  When  the  initial  size  of  the 
vortex  core  is  considerably  larger  than  the  sphere  size  (o'  ^  4),  the  effect  of  the  shear  flow 
(induced  by  the  passage  of  the  vortex  tube)  across  the  sphere  diminishes.  We  explained 
in  §3.2.1  that  the  magnitude  of  Q  depends  on  the  combined  effects  of  the 
downwash  and  the  shear  flow  across  the  sphere  due  to  the  vortex  tube.  As  a  result,  the 
magnitude  of  C,,.  ^  decreases  and  approaches  Q.  ,  as  o- 1>  1 .  and 

also  decrease  for  the  same  reason.  In  addition,  when  the  initial  size  of  the  vortex  core 
is  larger  than  the  sphere  size,  the  effect  of  the  wake  behind  the  sphere  on  the  vortex  tube 
diminishes.  As  a  consequence,  the  magnitudes  of  the  lift  and  moment  coefficients  decay 
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slowly  towards  zero  after  they  peak  near  the  time  of  passage  of  the  vortex  tube  centre 
by  the  plane  -  =  1. 

Summarizing  the  effects  of  the  vortex  size,  the  maximum  positive  lift  coefficient  and 
the  r.m.s.  lift  coefficient  depend  only  on  the  circulation  at  small  values  of  cr  while 
they  depend  only  on  (and  not  a)  at  large  values  of  cr.  For  mid-range  values  of  or, 
they  depend  on  both  <r  and  (or  equivalently  both  <r  and  F ^). 

In  §3.3,  we  investigated  the  effect  of  the  offset  distance  on  the  flow  field  for  Re  =  100 
and  o-  =  1 .  We  now  examine  the  effect  of  initial  offset  distance  of  the  vortex  tube  on 
the  lift  and  moment  coefficients  of  the  sphere  as  a  function  of  the  size  of  the  tube  at 
Re  =  100.  The  values  of  and  of  the  sphere  for  initial  offset 

distance  of  the  vortex  tube  in  the  range  of  -0.5 are  within  5% 
difference  from  their  values  for  d„ff  =  0.  On  the  other  hand,  ,nax  and  rm»  of  the 
sphere  for  initial  offset  distance  of  the  vortex  tube  in  the  same  range  vary  by  13  %  from 
their  values  for  d„fj  =  0. 


3.5.  Effects  of  Reynolds  number 

Computations  similar  to  those  in  §3.2  were  performed  for  four  different  Reynolds 
numbers  in  the  range  of  20  <  <  80,  d„ff  =  0,  and  o-  =  1  with  =  0.2  in  addition 

to  the  base  case  Re  =  100.  We  also  performed  the  same  calculation  with  two  different 
values  of  =  0.1  and  0.3  and  found  that  ^ L.maxV  ^L,max2^  ^L.rms’  ^M.max'  ^nd 
Cm.tm,  are  linearly  proportional  to  for  each  Reynolds  number.  Figures  (17H20) 
show  results  for  v^^x  —  0-2. 

Figure  17  shows  the  total  maximum  positive  lift  coefficient  and  the  coefficients  due 
to  pressure  and  viscous  contributions  as  a  function  of  Reynolds  number  for  o-  =  1  and 
=  0  with  v^ax  =  0.2.  The  coefficient  due  to  pressure  contribution  is  a  little  higher 
than  that  due  to  viscous  contribution  at  Reynolds  number  100.  Both  coefficients  due 
to  pressure  and  viscous  contributions  increase  as  the  Reynolds  number  decreases,  but 
the  viscous  coefficient  becomes  greater.  The  total  maximum  positive  lift  coefficient 
increases  with  a  negative  power  of  Reynolds  number  as  Reynolds  number  decreases, 
as  will  be  shown  in  figure  19  on  a  log-log  scale  and  is  expressed  by 

=  (23) 

where  the  constant  A  =  3.5  and  P  =  —0.32. 

Figure  18  shows  the  total  maximum  negative  lift  coefficient  and  the  coefficients  due 
to  pressure  and  viscous  contributions  as  a  function  of  Reynolds  number  for  o-  =  1  and 
=  0  with  =  0.2.  The  magnitude  of  the  coefficient  due  to  pressure  contribution 
is''2.38  times  higher  than  that  due  to  viscous  contribution  at  Reynolds  number  100.  As 
mentioned  in  §3.2,  the  sphere  experiences  the  maximum  negative  lift  coefficient  whose 
magnitude  is  greater  than  the  maximum  positive  lift  coefficient  when  the  vortex  tube 
passes  the  sphere,  because  the  vortex  tube  produces  high  fluid  velocity  gradient  across 
the  sphere  as  well  as  downwash  on  the  sphere.  (We  note  that  the  shear  flow  effect 
induced  by  the  vortex  tube  would  diminish  when  the  size  of  the  vortex  tube  becomes 
large.)  Thus,  the  pressure  contribution  is  much  higher  than  the  viscous  contribution  to 
the  total  maximum  negative  lift  coefficient.  The  magnitude  of  the  coefficient  due  to 
viscous  contributions  increases  as  the  Reynolds  number  decreases,  on  the  other  hand, 
that  due  to  pressure  contribution  decreases  as  the  Reynolds  number  decreases.  As  a 
consequence,  the  magnitude  of  the  total  maximum  negative  lift  coefficient  is  not 
sensitive  to  the  change  of  the  Reynolds  number  and  slowly  increases  as  the  Reynolds 
number  decreases. 
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Figure  17.  Total  maximum  positive  lift  coefficient  and  the  coefficients  due  to  pressure  and  viscous 
contributions  as  a  function  of  Reynolds  number  for  -  0  and  (r  =  1  with  =  0.2.  total;  0» 
pressure;  A»  viscous. 
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Figure  18.  Total  maximum  negative  lift  coefficient  and  the  coefficients  due  to  pressure  and  viscous 
contributions  as  a  function  of  Reynolds  number  under  the  same  conditions  as  figure  17.  □,  total;  O. 
pressure;  A.  viscous. 

Figure  19  shows  the  r.m.s.  lift  coefficient  and  also  the  maximum  positive  lift 
coefficient  as  a  function  of  Reynolds  number  on  a  log-log  scale  for  o-  =  1  and  =  0 
with  =  0.2.  The  r.m.s.  lift  coefficient  increases  with  a  negative  power  of  Reynolds 
number  as  the  Reynolds  number  decreases  and  is  expressed  by  (23)  with  A  =  2.3  and 
P  =  —0.275  for  cr  =  1.  The  effect  of  Reynolds  number  (20  ^  Re  ^  80)  on  the  lift 
coefficient  was  investigated  for  the  vortex  size  larger  than  o*  =  1  (2  ^  o'  ^  4).  The 
maximum  positive  lift  coefficient  and  the  r.m.s.  lift  coefficient  are  linearly  proportional 
only  to  and  independent  of  a  when  cr  ^  2  at  fixed  Reynolds  number  as  in  §3,4  for 
Re  =  100.  The  r.m.s.  lift  coefficient  is  expressed  by  (23)  with  A  =  SA  and  P  =  —0.45 
and  written  again  here  for  later  use. 

^L,  rms  ”  2  ^  (T  ^  4. 


(24) 
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Figure  19.  Root  mean  square  of  the  lift  coefficient  and  maximum  positive  lift  coefficient  as  a 
function  of  Reynolds  number  under  the  same  conditions  as  figure  17.  A,  maxi ;  O,  r.m.s. 
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Figure  20.  Maximum  moment  coefficient  and  root  mean  square  of  the  moment  coefficient  as  a 
function  of  Reynolds  number  under  the  same  conditions  as  figure  17.  A.  max;  O.  r.m.s. 


For  the  maximum  positive  lift  coefficient,  A  =  8.9  and  P  =  —0.45. 

Figure  20  shows  the  maximum  moment  coefficient  and  the  r.m.s.  moment  coefficient 
as  a  function  of  Reynolds  number  for  cr  =  1  and  d„ff  =  0  with  =  0.2.  The 
maximum  moment  coefficient  and  the  r.m.s.  moment  coefficient  are  affected  by  only 
the  viscous  effect  and  increases  with  a  negative  power  of  Reynolds  number  as  the 
Reynolds  number  decreases.  follow  the  form  of  (23)  with  ^  =  1 .95 

and  P  =  —0.56  for  the  former,  and  A  =  1.05  and  P  =  —0.665  for  the  latter  for  tr  =  1. 
The  effect  of  Reynolds  number  (20  ^  /?c  ^  80)  on  the  moment  coefficient  was 
investigated  for  2  ^  a-  <  4.  The  behaviour  of  the  moment  coefficient  at  each  Reynolds 
number  is  similar  to  that  of  the  moment  coefficient  at  Re  =  100  which  was  explained 
in  §3.4.  Cm.  max  and  follow  the  form  of  (23)  with  A  =  5.5  and  P  =  -0.83  for 

the  former,  and  A  =  3.1  and  P  =  —0.88  for  the  latter  for  2  ^  tr  ^  3. 
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Figure  21.  Lift  coefficients  of  the  sphere  as  a  function  of  time  and  for 
Re  =  20,  =  0  and  <r  =  1 . 


The  variation  of  the  lift  and  moment  coefficients  for  —0.5  ^  d„ff  <  from 

those  for  d^ff  =  0,  which  is  given  in  §3.4  for  Re  =  100,  decreases  at  fixed  cr  as  the 
Reynolds  number  decreases.  For  example,  at  Re  =  20,  the  difference  in  the  lift 
coefficient  is  4%  and  that  in  the  moment  coefficient  is  10%. 

We  investigate  the  lift,  moment  and  drag  coefficients  at  Reynolds  number  20  in  order 
to  find  out  about  the  Reynolds  number  effect  in  more  detail.  Figure  21  slows  the  lift 
coefficients  of  the  sphere  as  a  function  of  time  for  Re  =  20,  d„ff  =  0,  and  cr  =  1 .  The  lift 
coefficients  are  obtained  with  four  different  maximum  fluctuation  velocities  due  to  the 
vortex  tube,  =  01. 0-2,  0.3  and  0.4.  The  maximum  positive  lift  coefficient  Q 
occurs  at  t  =  6.6,  and  the  maximum  negative  lift  coefficient  occurs  at  t  =  12.5. 

The  lift  coefficient  is  linearly  proportional  to  the  maximum  fluctuation  velocity  (or  the 
circulation  of  the  vortex  tube)  at  each  time  over  the  whole  time  computed  (0  ^  r  ^  24). 
This  shows  that  the  nonlinear  effect  at  Re  =  20  is  much  less  than  that  at  Re  =  100.  In 
contrast  to  figure  13  which  shows  the  lift  coefficient  for  Re  =  100,  figure  21  shows  that 
the  lift  coefficient  decays  slowly  to  zero  after  it  attains  the  maximum  negative  value. 
This  indicates  that  for  Re  =  100  viscous  diffusion  in  the  wake  is  much  stronger  than 
that  in  the  upstream;  on  the  other  hand,  viscous  diffusion  is  uniformly  important  all 
around  the  sphere  at  Re  =  20  compared  to  Re  =  100.  The  behaviour  of  the  moment 
and  drag  coefficients  with  time  is  similar  to  that  of  the  case  of  Re  =  100.  The  time- 
averaged  lift  and  moment  coefficients  are  nearly  zero  (0(10"^)  and  0(10"®), 
respectively),  and  the  time-averaged  drag  coefficient  is  close  to  that  of  the  axisymmetric 
flow  without  the  vortex  tube  (the  difference  between  them  is  0(10"®)). 

One  of  the  reviewers  noted  that  some  of  our  results  can  be  explained  using 
dimensional  analysis  as  follows.  The  lift  force  on  the  sphere  can  be  expressed  in  a 
functional  form  as 

=AU'^,p\  a;,  r\  cr',  t’). 

From  dimensional  analysis,  it  follows  that 


Q(0  =  r 


f; 


=f\ 


U' 


r  'p'U'^a, 


02’  03’  04  ’ 
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where  are  the  four  parameters  appearing  on  the  left-hand  side  of  the  last  equal  sign. 
If  Reynolds  number  and  the  dimensionless  offset  distance  are  fixed,  Q(/)  will  be  a 
function  of  and  t.  Furthermore,  when  ar'/a'g  is  small,  we  expect  that  a'„  is  more 

important  than  <r'.  Then  should  be  a  function  only  of  4>i/4>2  and  t  yielding; 

In  the  opposite  limit,  a'^  should  be  unimportant  and 

Cflt)  =  =  AVmax ;  0-  (26) 

Equations  (25)  and  (26)  are  consistent  with  our  results  in  §3.4. 

Finally,  we  discuss  the  effect  of  the  initial  location  of  the  vortex  tube  upstream  from 
the  centre  of  the  sphere.  We  have  shown  earlier  that  the  maximum  positive  lift 
coefficient  is  expressed  at  given  cr  and  Re  as  „axi  —  ^^max  when  the  initial  location 
of  the  vortex,  /,  is  10  radii.  From  a  different  initial  location,  say  /*,  the  equation  for  the 
maximum  positive  lift  coefficient  should  be  modified  as  follows. 

^L,  max  l  ^^max 
_  O’*  ^ 


10.04(/*-/)Y'*  * 
Recr^  ) 


=  c*vt 


(27) 


where  we  used  =  <r*/<r  and  -<r  =  10.04  (/*  -1)1  Re,  which  are  obtained 

from  the  evolution  of  a  point  vortex  in  a  viscous  fluid  (Batchelor  1967).  denotes 
the  maximum  fluctuation  velocity  due  to  the  vortex  tube  whose  initial  location  is  /*  radii 
upstream  from  the  centre  of  the  sphere.  Note  that  the  proportionality  constant  c  is  now 
modified  as  c*  for  the  new  initial  location  of  the  vortex  /* 

We  examined  the  accuracy  of  (27)  by  performing  computations  for  /*  =  8  and  12 
with  the  same  parameters  as  used  in  §3.2. 1  except  the  initial  location  of  the  vortex  tube. 
The  magnitude  of  Q  obtained  from  (27)  differs  by  0.2%  from  that  of  the  full 
computations.  The  equation  of  the  r.m.s.  lift  coefficient  should  also  be  modified  as 
(27)  for  the  new  initial  location  of  the  vortex  tube.  In  addition,  the  time  span,  over 
which  averaging  the  lift  coefficient  is  performed  should  be  modified  according  to 

^ave  ‘at’«  '  v*  ^ 


4.  Conclusions 

As  a  first  step  towards  better  understanding  the  physics  of  interaction  between  a 
particle  and  the  turbulent  carrier  flow,  we  have  investigated  numerically  the  unsteady, 
three-dimensional,  incompressible,  viscous  flow  interactions  between  a  vortical 
(initially  cylindrical)  structure  advected  by  a  uniform  free  stream  and  a  spherical 
particle  suddenly  placed  and  held  fixed  in  space  for  a  range  of  particle  Reynolds 
numbers  20  ^  Re  ^  100.  The  counterclockwise  rotating  cylindrical  vortex  tube  is 
initially  located  ten  sphere-radii  upstream  from  the  centre  of  the  sphere. 
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A  summary  of  our  findings  and  their  applications  is  provided  as  follows. 

(i)  One  significant  finding  in  our  study  is  that  the  r.m.s.  lift  coefficient  for  a  particle 
is  linearly  proportional  to  the  upwash  (or  downwash)  induced  by  the  vortex  tube 
motion  on  the  particle  normal  to  the  direction  of  the  free  stream  in  our  case  (or  the 
direction  of  the  particle  motion  in  the  case  of  a  free  particle)  and  is  independent  of  the 
size  of  the  vortex  tube  when  the  size  of  the  vortex  is  greater  than  that  of  the  particle, 
2  <  O'  <  4.  This  result  can  be  applied  to  turbulent  flows  containing  small  concentration 
of  particles  in  order  to  obtain  the  r.m.s.  lift  force  on  a  particle.  A  turbulent  flow 
possesses  a  wide  spectrum  of  eddy  sizes.  The  large  eddies  contain  most  of  the  turbulent 
kinetic  energy  and  produce  high-velocity  fluctuations,  and  so  they  are  responsible  for 
the  dispersion  of  particles.  The  particle  size,  at  the  extremes,  may  be  comparable  to 
either  the  integral  lengthscale  or  to  the  Kolmogorov  lengthscale.  When  the  size  of 
particle  is  comparable  to  the  integral  lengthscale,  the  r.m.s.  lift  coefficient  of  the  particle 
is  obtained  by  (24).  Furthermore,  our  results  tend  to  support  the  idea  that  (24)  would 
be  applicable  to  the  case  of  an  eddy  much  larger  than  the  particle.  Thus,  when  the  size 
of  particle  is  comparable  to  the  Kolmogorov  lengthscale,  the  r.m.s.  lift  coefficient  of  the 
particle  can  be  calculated  approximately  by  (24),  where  **  ^he  maximum  velocity 
fluctuation  due  to  an  eddy  of  size  comparable  to  the  integral  lengthscale.  The  time 
during  which  the  particle  is  influenced  by  the  eddy  is  of  the  order  of  the  eddy  life  time. 

The  deflection  of  the  particle  path  will  depend  on  the  magnitude  of  the  r.m.s.  lift 
coefficient  and  the  ratio,  p,.,  of  the  particle  density  to  that  of  the  carrier  fluid 
(Q  =  where  A  is  the  dimensionless  acceleration  of  the  particle  due  to  the  lift 
force).  This  result  provides  a  simple  method  to  estimate  the  deflection  of  particle 
trajectory  in  the  dilute  particle-laden  turbulent  flow.  Equation  (24)  and  the  non- 
dimensionalized  Newton’s  second  law  show  that  the  deflection  increases  slowly  as  the 
Reynolds  number  of  the  particle  decreases. 

(ii)  The  magnitude  of  the  r.m.s.  moment  coefficient  of  the  particle  is  one  order  of 
magnitude  less  than  that  of  the  r.m.s.  lift  coefficient  when  Re  ^  20.  Furthermore,  when 
the  initial  size  of  the  vortex  core  is  considerably  larger  than  the  sphere  size  (cr  ^  4),  the 
effect  of  the  shear  flow  (induced  by  the  passage  of  the  vortex  tube)  across  the  sphere 
diminishes  and  the  torque  on  the  particle  decreases.  Thus,  the  torque  on  the  particle 
might  be  negligible  in  many  applications. 

(iii)  When  a  vortex  tube  advected  by  a  uniform  free  stream  approaches  a  sphere,  the 
sphere  experiences  lower  drag  than  that  of  a  sphere  subjected  to  an  axisymmetric  flow 
if  the  sphere  were  free  to  move  rather  than  fixed,  unless  the  initial  offset  distance  of  the 
vortex  tube  is  large  positive,  as  explained  in  §3.3.2.  The  lower  drag  is  caused  by  the 
upward  motion  of  the  sphere  due  to  the  upwash  of  the  approaching  vortex  tube,  and 
thus  the  centre  of  the  vortex  tube  would  be  located  below  the  front  stagnation  point 
of  the  sphere.  This  causes  lower  dynamic  pressure  ahead  of  the  front  stagnation  point. 

(iv)  Some  interesting  unsteady  phenomena  in  the  near  wake  have  been  discovered. 
The  shape  of  the  near  wake  behind  the  spherical  particle  is  controlled  by  the  pressure 
difference  between  the  top  and  bottom  of  the  near  wake  as  was  indicated  by  Kim  et  al. 
(1993).  The  instantaneous  flow  patterns  around  a  spherical  particle  in  a  turbulent  flow 
would  include  some  of  those  described  in  this  paper.  For  example,  our  recent  results 
(to  be  published),  from  a  study  of  the  interactions  between  a  vortex  pair  advected  by 
a  uniform  free  stream  and  a  sphere,  show  that  the  streamlines  are  similar  to  those 
described  in  the  present  paper. 

We  are  also  studying  the  heat  and  mass  transfer  for  a  droplet  interacting  with  an 
array  of  vortices  in  high-temperature  and  high-pressure  environments  such  as  that  in 
a  gas  turbine  combustor. 
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Abstract 

This  study  concerns  the  detailed  three-dimensional, 
viscous,  incompressible  interactions  of  vortical  struc¬ 
tures  with  a  rigid  sphere.  A  pair  of  vortical  structures 
(initially  cylindrical)  advect  past  a  sphere.  Navier- 
Stokes  equations  describe  the  unsteady  viscous  flow 
field.  Finite-difference  computations  yield  flow  prop¬ 
erties  plus  temporal  behavior  of  lift,  drag,  cind  moment 
coefficients  of  the  sphere.  Lift  and  moment  coefficients 
are  shown  to  be  linearly  proportional  to  the  maximum 
velocity  fluctuation.  Effects  of  Reynolds  number,  initial 
vortex  size,  and  initial  vortex  configuration  are  deter¬ 
mined.  Lift  coefficients  are  used  to  estimate  spherical 
particle  deflection  in  turbulent  flows;  deflection  is  found 
to  increase  slowly  as  Reynolds  number  decreases.  The 
moment  coefficient  is  an  order  of  magnitude  less  tham 
the  lift  coefficient  implying  that  torque  is  often  negligi¬ 
ble. 


Nomenclature 


^£f,maxl 

CL,max2 

CL^rm$ 

doff 

Ni,N2,Ns 

Re 

u,  u,  u; 

^max 

‘^maxi 

UL 

t 


dimensional  sphere  radius 
maximum  positive  lift  coefficient 
maximum  negative  lift  coefficient 
rms  lift  coefficient 

initial  offset  distance  normalized  by  a' 
numbers  of  grids  in  r),  ^  directions 
Reynolds  number  based  on  sphere  diameter 
velocities  in  the  x,  y,  and  z  direction 
maximum  fluctuation  velocity  induced 
by  one  vortex  tube 

total  maximum  fluctuation  velocity  induced 
by  a  pair  of  vortex  tubes 
dimensional  freestream  velocity 
time  normalized  by 
Cartesian  coordinates 
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Greek  symbols 

Donorthogonal  generalized  coordinates 
(T  radius  of  vortex  tube  normalized  by 

circulation  of  vortex  tube  normalized 
by  and  aj, 

y-component  of  vorticity  vector 

Superscript 

^  dimensional  quantity 

Subscript 

o  initial  quantity 

1.  Introduction 

Particle-laden  turbulent  flows  are  important  in  many 
natural  and  engineering  applications  such  as  atmo¬ 
spheric  dispersion  of  pollutants,  combustion  systems, 
and  transport  of  suspended  substances  in  slurries  or 
pneumatic  systems.  The  fluid  may  be  either  gas  or 
liquid  containing  particles,  droplets,  or  bubbles.  The 
geometric  scale  of  the  flow  may  be  as  large  as  weather 
patterns  in  the  atmosphere  of  a  planet,  or  as  small  as 
the  veins  transporting  blood  ceils  in  the  body’s  circu¬ 
latory  system. 

In  general,  there  is  no  analytical  solution  to  the  equar 
tion  of  motion  of  a  single  particle  in  a  laminar  or  tur¬ 
bulent  flow  due  to  the  nonlinear  nature  of  the  equation. 
This  has  led  to  extensive  application  of  numerical  ap¬ 
proaches  to  study  the  particle  motion  in  a  turbulent 
flow.  When  the  concentration  of  particles  in  a  flow  is 
small,  the  particle-peirticle  interaction  such  as  collisions 
or  repulsions  can  be  neglected,  and  the  flow-particle 
interaction  is  dominated  by  the  time  response  of  the 
individual  particles.  This  type  of  flow  is  classified  as 
a  dilute  dispersed  flow.  Extensive  research  has  been 
performed  to  characterize  a  variety  of  dilute  dispersed 
flow  fields.  Review  articles  by  Lumley  [1],  Crowe  [2], 
and  Faeth  [3]  provide  a  summary  of  current  research 
results  and  directions. 

Most  studies  to  predict  particle  motion  in  a  tur¬ 
bulent  flow  employ  the  Eulerian-Lagrangian  approach, 
where  particles  are  considered  as  point  sources  in  the 
Navier-Stokes  equations.  After  the  local  fluid  veloc¬ 
ities  are  obtained  from  the  Navier-Stokes  equations, 
the  particle  motion  is  obtained  by  solving  the  BBO 
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(Basset-Boussinesq-Oseen)  equation.  The  drag  force  in 
that  equation  consists  of  two  terms,  namely,  the  quasi¬ 
steady  Stokes  drag  and  the  Basset  correction.  Strictly, 
this  correction  is  only  correct  analytically  as  a  correc¬ 
tion  to  the  linear  theory  developed  for  Stokes  flow  and 
valid  only  for  very  low  particle  Reynolds  number.  Thus, 
it  is  not  possible  in  that  approach  to  obtain  the  actual 
force  2uid  torque  on  the  particle,  heat  and  mass  transfer 
rate  of  the  particle,  and  the  details  of  the  flow  around 
the  particle  including  near  wake  due  to  the  turbulent 
velocity  fluctuation  in  a  finite-Reynolds-number-flow. 
An  alternative  is  to  examine  the  details  of  the  inter¬ 
action  between  a  particle  and  its  surrounding  turbu¬ 
lent  flow  by  solving  numerically  the  full  Navier-Stokes 
equations  around  an  individual  particle  with  the  ap¬ 
propriate  boundary  conditions  and  employing  a  simple 
mathematical  description  for  turbulent  velocity  fluctu¬ 
ation. 

A  turbulent  flow  possesses  a  wide  spectrum  of  eddy 
sizes.  For  example,  in  a  typical  gas  turbine  combus¬ 
tor  where  the  Reynolds  number  is  of  the  order  of  10^ 
and  the  integral  length  scale  is  of  the  order  of  0.1m, 
the  smallest  (Kolmogorov)  length  scale,  is  about 
lOO/im,  which  is  comparable  to  the  size  of  a  typical  fuel 
droplet.  Motion  at  the  large  (integral)  length  scale  con¬ 
tains  most  of  the  turbulent  kinetic  energy  and  governs 
the  dispersion  of  particles  (or  droplets).  On  the  other 
hand,  fluid  motion  at  the  Kolmogorov  length  scale  ex¬ 
periences  the  largest  strain  rates  and  scalar  gradients 
in  the  flow.  The  largest  scalar  gradients  control  the 
important  phenomena  of  heat  and  mass  transfer  and 
chemical  reaction. 

As  a  first  step  towards  better  understanding  the  in¬ 
teractions  of  a  particle  with  a  turbulent  flow,  Kim,  El- 
ghobashi  ic  Sirignano  [4]  studied  the  three-dimensional 
flow  interactions  between  a  vortical  (initially  cylindri¬ 
cal)  structure  advected  by  the  free  stream  and  a  spher¬ 
ical  solid  particle  held  fixed  in  space.  The  particle 
Reynolds  number  based  on  the  freestream  velocity  and 
the  particle  diameter  was  in  the  range  20  <  iZe  <  100. 
The  initial  size  of  the  cylindrical  vortex  tube  was  in  the 
range  0.25  <  <  4,  where  er  is  the  radius  of  the  vortex 

tube  normalized  by  that  of  the  particle.  They  found 
that  the  rms  lift  coefficient  of  the  sphere  is  linearly  pro¬ 
portional  to  the  circnlaiion  of  the  vortex  tube  at  small 
values  of  o’,  on  the  other  hand,  at  large  values  of  <r,  the 
rms  lift  coefficient  is  linearly  proportional  to  the  max¬ 
imum  fluctuation  velocity  induced  by  the  vortex  tube 
but  independent  of  a. 

We  study  in  this  paper  the  interactions  between  a 
pair  of  vortex  tubes  and  a  rigid  sphere.  Our  main  goal 
is  to  generalize  the  findings  of  Kim  ei  al  [4]  by  deter¬ 
mining: 

1.  the  details  of  the  interaction  of  a  pair  of  vortex 
tubes  with  each  other  and  with  the  sphere  in  a 
viscous  incompressible  flow 

2.  the  relationship  between  the  lift  coefficient  of  the 


sphere  and  the  maximum  fluctuation  velocity  in¬ 
duced  by  a  pair  of  vortex  tubes 

3.  the  effects  of  Reynolds  number,  vortex  size,  and 
initial  offset  distance  of  the  vortex. 

The  detailed  study  of  the  interactions  between  the 
particle  and  the  unsteady  velocity  field  provides  funda¬ 
mental  information  about  the  flow  behavior  which  can 
be  used  in  developing  mathematical  models  for  two- 
phase  flows. 

2.  Formulation  and  numerical  solution 
2.1  The  flow  description 

Consider  the  time-dependent,  three-dimensional,  in¬ 
compressible,  viscous  flow  interactions  between  a  pair 
of  initially  cylindrical  vortex  tubes  and  a  solid  sphere. 
The  vortex  tubes  are  moving  with  the  launinar  free 
stream,  and  a  sphere  is  suddenly  placed  and  held  fixed 
in  space  as  shown  in  Figure  1.  The  initial  offset  dis¬ 
tance,  <//.  denotes  the  shortest  distance  between  the 
initial  vortical  axis  and  the  y-z  plane  which  is  parallel 
to  the  free  stream.  All  the  variables  are  nondimension- 
alized  using  the  sphere  radius  as  the  characteristic 
length  and  [7^  as  the  characteristic  velocity,  where  the 
superscript  /  denotes  dimensional  quantity.  The  pair 
of  cylindrical  vortex  tubes,  whose  diameters  are  of  the 
same  size  and  of  the  order  of  the  sphere  diameter,  are 
initially  located  ten  radii  upstream  from  the  center  of 
the  sphere.  The  effects  of  the  vortex  tube  on  the  sphere 
are  negligible  at  this  initiad  distance  because  the  mag¬ 
nitude  of  the  initial  velocity  field  induced  by  the  vortex 
tubes  is  less  than  2  percent  of  the  free  stream  veloc¬ 
ity.  Far  upstream,  the  flow  is  uniform  with  constant 
velocity  parallel  to  the  y-z  plane.  We  have  one 
symmetry  plane,  the  x-z  plane,  as  seen  in  Figure  1.  A 
second  symmetry  plane  (y-z)  exists  when  the  two  vor¬ 
tices  have  opposite  rotations.  Our  generzJ  formulation 
does  not  take  advantage  of  this  second  symmetry. 

Two  coordinate  systems  are  used  in  our  formulation: 
the  Cartesian  coordinates  (x,y,z)  and  the  nonorthogo 
nal  generalized  coordinates  (^,j;,C)-  The  origin  of  the 
former  coincides  with  the  sphere  center.  (  is  the  radial, 
i;  is  the  angular,  and  (  is  the  azimuthal  coordinates. 
The  nonorthogonai  generalized  coordinate  system  can 
be  easily  stdapted  to  three-dimensional  arbitrary  ge¬ 
ometries.  In  the  present  study,  a  spherical  domain  is 
used,  and  the  grid  reduces  to  an  orthogonal,  spheri¬ 
cal  grid.  The  grids  are  denser  near  the  surface  of  the 
spherical  particle,  and  the  grid  density  in  the  radial 
direction  is  controlled  by  the  stretching  function  devel¬ 
oped  by  Vinokur  [5].  Due  to  symmetry,  the  physical 
domain  is  reduced  to  a  hadf  spherical  space.  The  do¬ 
main  of  the  flow  is  bounded  by  1  <  ^  <  ,  1  <  ry  < 

1  <  C  <  Wa,  where  {  =  1  and  Ni  correspond,  respec¬ 
tively,  to  the  sphere  surface  and  the  farfield  boundary 
surrounding  the  sphere;  ?;  =  1  and  denote,  respec¬ 
tively,  the  positive  z-axis  (upstream)  and  the  negative 
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z-axis  (downstream);  C  =  1  refer,  respectively, 

to  the  x-z  plane  in  the  positive  x-direction  and  the  x- 
z  plane  in  the  negative  x-direction.  Uniform  spacing 
=  5?;  =  5^  =  1)  is  used,  for  convenience,  for  the 
generalized  coordinates. 

The  initial  vortex  tube  has  a  small  core  region  with 
a  radius  c  (normalized  by  the  sphere  radius).  This  core 
is  defined  such  that  the  initial  velocity  induced  by  the 
vortex  tube  approaches  zero  as  the  distance  from  the 
center  of  the  vortex  tube  goes  to  zero,  and  at  distances 
much  greater  than  tr,  the  induced  velocity  approaches 
that  of  a  point  vortex.  We  use  the  vortex  tube  con¬ 
struction  of  Spalart  [6],  which  has  the  following  stream 
function: 

Z,  i  =  0)  =  -^/n[(i  -  Xjf  +  (z-  Zj)^  +  a^] , 

where  Fy  is  the  nondimensional  circulation  around  the 
vortex  tube  at  radius  <r  and  at  the  initial  time.  Fj  is 
positive  when  the  vortex  tube  rotates  counterclockwise, 
and  Xj  and  Zj  denote  the  location  of  the  center  of  the 
vortex  tube.  The  circulation  around  a  circular  path  far 
away  from  the  center  of  the  vortex  is  given  by  Tjt  = 
2Fj.  The  stream  function  for  a  pair  of  vortex  tubes  is 
given  by 

z,  f  =  0)  =  -  -  xj)^  +  (z  -  z,- 

i=i 

(1) 

2.2  Governing  equations  and  boundary  condi¬ 
tion 


The  continuity  and  momentum  equations  and  the  ini¬ 
tial  and  boundary  conditions  are  nondimensionalized 
using  the  sphere  radius  as  the  characteristic  length 
and  as  the  characteristic  velocity. 

V-K  =  0  (2) 

^  +  wv  =  -Vp  +  w .  (3) 


The  governing  equations  (2)  and  (3)  are  cast  in  terms 
of  the  generalized  coordinates  (^ ,  j/,  Q  to  treat  a  three- 
dimensional  body  of  arbitrary  shape.  The  numerical 
integration  is  performed  using  a  cubic  computational 
mesh  with  equal  spacing  (i^  =  tfr;  =  =  !)• 

The  velocities  on  the  sphere  surface  are  zero  due  to 
the  no-slip  condition,  and  the  pressure  on  the  sphere  is 
obtained  from  the  momentum  equation.  The  boundary 
conditions  are 


dp  _  2  d^Vn 
dn  Re  dri^  ’ 


u  =  v  =  ti;  =  0at^=rl, 


(4) 


p=0,  a  =  w  =  0,  w=la,t^  =  Ni,  N2m  <  »?  <  ^2 

(5) 


P  =  0.  ^  =  ^  =  ^=  0atf  =  iVi,l<»?<  N2m 

(6) 

dp  du  dw 


where  «,  v,  and  w  are  the  velocities  in  the  x,  y,  and  z 
direction,  respectively,  is  the  velocity  in  the  direc¬ 
tion  normal  to  the  sphere  surface,  and  p  is  the  pressure, 
n  denotes  the  direction  normal  to  the  sphere  surface, 
d/dn  =:  and  rj  =  N^m  denotes  the 

mid-plane  between  17  =  1  and  W2,  Equation  (7)  corre¬ 
sponds  to  the  symmetry  conditions  and  zero  v  velocity 
in  the  x-z  symmetry  plane.  Conditions  guaranteeing 
continuity  in  the  q  direction  for  p,  u,  w,  and  w  on  the 
axes  17  =  1  and  17  =  W2  are  also  imposed. 

In  order  to  start  the  numerical  solution  of  equations 
(2)  and  (3),  we  provide  initial  velocity  by  superposing 
the  flow  fields  due  to  the  uniform  stream  and  the  vortex 
tube  and  the  no-slip  condition  on  the  sphere  surface: 


Po  =  0,  Uo  = 


dz  ’ 


«.  =  0,  =  1  +  ^ 
except  at  ^  =  1 


(8) 


Po  =  0,  =  Vo  =  u^o  =  0  at  f  =  1,  (9) 

where  V^^  is  given  by  equation  (1). 

The  only  nondimensional  groupings  appearing  in  the 
equations  and  initial  and  boundary  constraints  are  the 
sphere  Reynolds  number,  vortex  tube  radius,  offset  dis¬ 
tance,  and  vortex  circulation  (or  vortex  Reynolds  num¬ 
ber). 

The  drag,  lift,  and  moment  coefficients  are  evaluated 
in  dimensional  form  as  follows. 


F'o  =  J^-p'n-k  dS'  +  J  n-r'-kdS'  (10) 


—p^n*i  dS*  +  J  dS*  (11) 
M'  =  J  t'  XT'  dS',  (12) 

where  5^  denotes  the  surface  of  the  sphere,  n  is  the 
outward  unit  normal  vector  at  the  surface,  r'  is  the  po¬ 
sition  vector  from  the  center  of  the  sphere,  and  r*  is  the 
viscous  stress  tensor.  The  lift  force  is  assumed  positive 
when  it  is  directed  toward  the  positive  x-axis.  Due  to 
symmetry,  only  the  y-component  of  the  moment  is  non¬ 
zero  and  is  assumed  positive  in  the  counter-clockwise 
direction. 

The  nondimensional  coefficients  of  drag,  lift,  and  mo¬ 
ment  are  defined  r^pectively  as 
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Cd  = 


Pd 


(13) 


3*2  Interactionfl  of  a  pair  of  vortex  tubes  of  like 
rotation  and  a  sphere 


Cl  = 


n 


(14) 


Cm  — 


M'^j 


2.3  Numerical  solution 


(15) 


A  three-dimensional,  implicit,  finite-difference  algo¬ 
rithm  has  been  developed  to  solve  simultaneously  the 
set  of  the  discretized  partial  differential  equations.  The 
method  is  based  on  an  Alternating-Direction-Predictor- 
Corrector  (ADPC)  scheme  to  solve  the  time-dependent 
Navier-Stokes  equations.  ADPC  is  a  slight  variation  of 
Alternating-Direction-Implicit  (ADI)  method  and  im¬ 
plemented  easily  when  embedded  in  a  large  iteration 
scheme  (Patnaik  [7]).  The  control  volume  formulation 
is  used  to  develop  the  finite-difference  equations  from 
the  governing  equations  with  respect  to  the  general¬ 
ized  coordinates  An  important  part  of  solving 

the  Navier-Stokes  equations  in  primitive  variables  is  the 
calculation  of  the  pressure  field.  In  the  present  work, 
a  pressure  correction  equation  is  employed  to  satisfy 
indirectly  the  continuity  equation.  The  pressure  cor¬ 
rection  equation  is  of  the  Poisson  type  and  is  solved  by 
the  Successive-Over-Relaxation  (SOR)  method. 

The  overall  solution  procedure  is  based  on  a  cyclic 
series  of  guess-sind-correct  operations.  The  velocity 
components  are  first  calculated  from  the  momentum 
equations  using  the  ADPC  method,  where  the  pressure 
field  at  the  previous  time  step  is  employed.  This  esti¬ 
mate  improves  as  the  overall  iteration  continues.  The 
pressure  correction  is  calculated  from  the  pressure  cor¬ 
rection  equation  using  the  SOR  method,  and  new  es¬ 
timates  for  pressure  and  velocities  are  obtained.  This 
process  continues  imtil  the  solution  converges  at  each 
time  step. 


3.  Results  and  discussion 


3.1  Numerical  accuracy 


Here  we  examine  the  flow  generated  by  an  impul¬ 
sively  started  solid  sphere  in  a  quiescent  fluid  at  two 
Reynolds  numbers:  20  and  100.  The  time-dependent 
solution  converges  asymptotically  to  a  steady-state 
which  is  in  good  agreement  with  the  available  exper¬ 
imental  data  and  correlations.  Refer  the  details  for  nu¬ 
merical  data  to  Kim  ei  ai  [4].  Although  the  solution 
in  these  test  cases  are  axis3mimetric,  none  of  the  three 
velocity  components  in  our  formulation  becomes  iden¬ 
tically  zero.  Therefore,  the  three-dimensional  solution 
scheme  is  fully  exercised  here.  The  51  x  51  x  51  grid  is 
used  in  the  following  calculations. 


Now  consider  in  particular  the  interactions  of  a  pair 
of  vortex  tubes  advected  by  the  free  stream  and  a 
sphere  suddenly  placed  in  the  flow  and  held  fixed  in 
space.  The  two  vortex  tubes  are  initially  of  the  same 
size  and  cylindrical  shape  rotating  counter-clockwise  as 
shown  in  Figure  1.  The  two  vortex  tubes  are  separated 
by  the  same  distance  from  the  y-z  plane  so  that  the 
offset  distance  of  one  vortex  tube  is  minus  the  offset 
distance  of  the  other.  The  center  of  the  each  vortex 
tube  is  located  at  10  sphere-radii  upstream  from  x-y 
plane  containing  the  center  of  the  sphere.  The  base 
case  calculation  is  that  of  Re  =  100,  do//  =  ±1.5,  and 

(T  =  1. 

Initially  each  vortex  tube  has  the  maximum  fluctu¬ 
ation  velocity  v^ax  on  the  edge  of  the  core.  Because 
the  velocity  and  vorticity  fields  induced  by  one  vor¬ 
tex  tube  influence  those  by  the  other,  the  total  maxi¬ 
mum  fluctuation  velocity,  Vmaxti  induced  by  two  vor¬ 
tex  tubes  depends  on  the  separation  distance  between 
them  and  their  size  and  is  in  the  range  Vmax  <  ^maxt  ^ 
2vma*  with  2  Vmax  when  |do//|  =  0  and  Vmax  when 
|do//|  »  1.  For  example,  Vmaxt  is  0.59  for  Vmax  =  0.4, 
do//  =  ±1.5,  and  <r  =  1. 

The  x-z  plane  of  symmetry  is  defined  as  the  principal 
plant,  where  the  strongest  interactions  occur  between 
the  vortical  structure  and  the  sphere.  The  line  connect¬ 
ing  the  front  and  rear  stagnation  points  in  the  standard 
axisymmetric  flow  over  a  single  sphere,  which  is  the  x 
=  0  line  in  the  principal  plane,  will  be  used  as  a  refer¬ 
ence  line.  We  refer  to  the  region  above  the  line  as  the 
‘upper’  region  and  that  below  the  line  as  the  ‘lower’ 
region. 

Figures  2(a)-(f)  display  the  contour  lines  of  y- 
component  vorticity  in  the  principaJ  plane  at  t  =  1, 
6,  10,  15,  21,  and  30  for  Re  =  100,  d^jf  =  ±1.5,  cr  =  1 
with  Vmaxt  =  0.59.  The  contour  values  of  the  vorticity 
are  ±0.4,  ±0.8,  ±1.4,  ±2,  with  the  highest  magnitude  at 
the  sphere  surface.  The  solid  and  dotted  lines  in  the 
figures  represent  positive  and  negative  values,  respec¬ 
tively. 

The  vorticity  contours  in  figures  2(a)  and  2(b)  show 
that  the  vortex  tubes  not  only  are  advected  down¬ 
stream  but  also  rotate  each  other.  The  contour  lines  of 
vorticity  in  the  figures  also  show  that  the  viscous  diffu¬ 
sion  takes  place.  Figures  2(c)  and  2(d)  show  that  the 
vortex  tubes  contact  the  boundary  layer  of  the  sphere 
and  go  around  the  bottom  of  the  sphere.  When  the 
vortex  tube  comes  close  to  the  sphere  boundary  layer, 
it  augments  the  magnitude  of  the  edge  velocity  in  the 
lower  boundary  layer  and  reduces  the  edge  velocity 
in  the  upper  boundary  layer.  The  result  is  a  higher 
strength  vorticity  in  the  lower  boundary  layer  than  in 
the  upper  boundary  layer  (see  the  vorticity  contours  in 
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figure  2(c)).  Consequently,  the  vorticity  in  the  lower 
boundary  layer  induces  a  velocity  in  the  downward  di¬ 
rection  at  the  location  of  the  vortex  tube  with  higher 
magnitude  than  that  induced  by  the  vorticity  in  the  up¬ 
per  boundary  layer.  This  downward  induced  velocity 
advecta  the  vortex  tube  below  the  sphere  (Kim  ei  oL 
[4]).  Figure  2(e)  shows  that  the  pairing  vortex  tubes 
merge  into  one  vortex  due  to  the  attraction  of  each 
other  and  the  viscosity. 

3.2.1  Lift,  moment^  and  drag  coefficients  and 
effect  of  tube  circulation 

Figure  3  shows  the  lift  coefficients  of  the  sphere  as  a 
function  of  time  for  Re  =  100,  =  ±1.5,  and  <r  =  1. 

The  lift  coefficients  are  computed  for  four  different  total 
maximum  fluctuation  velocities  v^axt  induced  by  the 
pair  of  vortex  tubes,  with  magnitudes  equal  to  0.148, 
0.295,  0.443,  and  0.590  (w  =  0.1,  0.2,  0.3,  and  0.4). 
The  lift  coefficient  of  the  sphere  interacting  with  a  sin¬ 
gle  vortex  tube  as  a  function  of  time  is  also  shown  as  a 
reference  in  figure  3  for  iZe  =  100,  do//  =  0,  and  cr  =  1 
with  Vjnaxi  =  0.148.  The  connection  between  the  case 
of  a  single  vortex  tube  and  that  of  a  pair  of  vortex  tubes 
will  be  discussed  in  the  next  section  in  detail.  Due  to 
the  sudden  placement  of  the  sphere  into  the  stream,  it 
takes  a  short  time  (0  <  f  <  0.8)  for  the  initial  flow 
perturbations  to  vanish. 

When  the  pair  of  vortex  tubes  approach  the  sphere 
(0  <  t  <  9),  they  produce  upwash  resulting  in  a  pos¬ 
itive  lift  force  on  the  sphere.  The  maximum  positive 
lift  coefficient  Cl, maxi  occurs  at  about  t  =  6.8.  On 
the  other  hand,  when  the  vortex  tubes  pass  the  sphere, 
they  produce  downwash  and  higher  fluid  velocity  near 
the  bottom  of  the  sphere  than  the  top  due  to  the  shear 
flow  imposed  by  the  vortex  tubes  resulting  in  a  negative 
lift  force.  The  magnitude  of  the  negative  lift  is  greater 
than  the  positive  lift.  The  maximum  negative  lift  coef¬ 
ficient  Climax!  occurs  at  about  t  =  12.2.  Ci^maxi  and 
CL,max2  are  linearly  proportional  to  the  maximam  fluc¬ 
tuation  velocity.  The  maximum  positive  lift  coefficient 
Cl, maxi  ifi  expressed  by 

Cl, maxi  ~  ^  Vnioj?!  ,  (1^) 

where  the  proportionality  constant  c  =  0.88.  The  max¬ 
imum  negative  lift  coefficient  CL,max2  is  also  expressed 
by  equation  (16)  but  with  c  =  —1.62.  After  the  lift  co¬ 
efficient  reaches  its  maximum  negative  value,  it  decays 
quickly  towards  zero  because  the  vortex  tube  vorticity 
is  diffused  in  the  sphere  wake.  The  time  averaged  lift 
coefficient  (averaged  over  a  time  span  between  i  =  0.8 
and  the  maximum  time  24.5)  for  all  values  of  Vmaxt  i® 
nearly  zero  (0(10“^)).  As  mentioned  earlier,  the  be¬ 
havior  of  CL{i)  during  the  period  0  <  t  <  0.8  is  in¬ 
fluenced  by  the  initial  flow  perturbation,  and  thus  its 
value  during  this  initial  period  is  excluded  from  the  av¬ 
eraging  process.  The  root  mean  square  CL,rm»  of  the 


lift  coefficient  as  a  function  of  time  is  also  linearly  pro¬ 
portional  to  Vjnaxt  C  =  0.7. 

Figure  4  shows  the  temporal  development  of  the  mo¬ 
ment  coefficients  for  the  sphere  under  the  same  con¬ 
ditions  of  Figure  3.  The  moment  coefficient  of  the 
sphere  interacting  with  a  single  vortex  tube  as  a  func¬ 
tion  of  time  is  also  shown  as  a  reference  in  figure  4  for 
Re  =  100,  dojj  =  0,  and  cr  =  1  with  Vmaxt  =  0.148. 

When  the  vortex  tubes  pass  the  sphere,  the  front 
stapation  point  on  the  sphere  is  shifted  above  the 
plane  a:  =:  0  due  to  the  downwash.  This  causes 
higher  shear  stress  in  the  lower  left  region  compared  to 
the  upper  left  region  resulting  in  a  positive  (counter¬ 
clockwise)  torque  on  the  sphere.  The  upward  shift  of 
the  ftont  stagnation  point  also  causes  the  shear  stress 
to  be  higher  in  the  top  and  upper  right  regions  than  in 
the  bottom  and  lower  right  regions  resulting  in  a  neg¬ 
ative  torque  on  the  sphere.  However,  the  effect  of  this 
negative  torque  is  diminished  by  the  shear  flow  induced 
by  the  vortex  tubes  across  the  sphere  which  produces 
high  shear  stress  at  the  bottom  of  the  sphere.  As  a  con¬ 
sequence,  a  net  high  positive  torque  acts  on  the  sphere. 
The  maximum  positive  moment  coefficient  Cm, mar  oc¬ 
curs  at  f  =:  11.5.  Cm, max  is  linearly  proportional  to 
^maxt  with  C  =  0,11. 

When  the  vortex  tubes  approach  the  sphere  or  aire 
relatively  far  away  from  the  sphere,  the  effect  of  the 
shear  flow  induced  by  the  vortex  tubes  across  the  sphere 
is  small,  resulting  in  a  net  weak  torque  on  the  sphere. 

The  time  averaged  moment  coefficient  (averaged  over 
a  time  span  between  t  =  0.6  and  24.5)  for  all  values  of 
Vmaxt  0(10“^),  and  the  root  mean  square  CM,rms  of 
the  moment  coefficient  is  linearly  proportional  to  Vmaxt 
with  c  =  0.043.  We  note  that  the  torque  depends  only 
the  distribution  of  the  shear  stresses  (r,.#  and  r^^)  and 
is  relatively  small  compared  to  the  lift  force. 

Figure  5  shows  the  drag  coefficients  of  the  sphere  as 
a  function  of  time  for  the  same  conditions  of  Figure 
4.  The  drag  coefficients  are  computed  for  four  different 
values  of  as  in  Figure  4,  in  addition  to  Vmaxt  =  0 
which  corresponds  to  the  axisymmetric  flow  without 
the  vortex  tubes. 

The  initially  large  values  of  shear  stress  and  pre^ure 
on  the  sphere  results  in  large  drag  as  shown  in  Figure  5. 
When  the  vortex  tubes  approach  the  sphere,  the  pres¬ 
sure  at  the  front  stagnation  point  is  lower  than  that 
of  the  axisymmetric  flow  past  a  sphere  due  to  the  low 
pressure  at  the  center  of  the  vortex  tube.  Also,  the 
maximum  shear  stresses  in  the  upper  and  lower  regions 
of  the  sphere  are  lower  than  those  of  the  axisymmetric 
flow.  This  causes  the  drag  on  the  sphere  to  be  lower 
than  that  of  the  axisymmetric  flow  without  the  vortex 
tube.  As  the  vortex  tubes  move  airound  the  bottom  of 
the  sphere,  the  front  stagnation  point  is  shifted  above 
the  plane  a?  =  0  due  to  the  downwash.  Consequently, 
high  pressure  and  high  shear  stress  act  respectively  in 
the  upper  left  region  and  the  lower  left  region.  This 


5 


increases  the  drag  during  the  period  9  <i  <  13.4.  For 
t  >  13.4,  the  drag  approaches  that  of  the  axisymmet- 
ric  flow  as  the  vortex  tube  moves  further  downstream. 
The  time  averaged  value  of  the  deviation  of  the  drag 
coefficient  from  that  of  the  axisymmetric  flow  past  a 
sphere  for  all  values  of  Umo*  is  nearly  zero  (0(10“''*)). 

3.2.2  Effects  of  the  size  and  the  offset  diatanre 
of  the  vortex  tube 

The  effects  of  the  size  of  the  vortex  tube  on  the 
flow  field  are  studied  by  performing  computations  sim¬ 
ilar  to  those  in  section  3.2.1  for  Re  =  100,  dp//  =  ±1.5, 
and  five  different  sizes  of  the  vortex  tube,  0.25  <0"  <4 
in  addition  to  the  base  case  (T  =  1. 

Table  1  shows  Ct.marif  Cl  ,max3t  Cl  ^rmsf  Cm  .mox) 
^d  CMrm*  ^  ^  function  of  the  vortex  tube  size  which 
covers  six  different  initial  radii  of  the  vortex  tube, 

(T  =  4,3,2, 1,0.5,  and  0.25,  for  Vmaxt  =  0.1.  Another 
computation  with  different  v^axt  showed  that  all  the 
lift  and  moment  coefficients  are  linearly  proportional  to 
Vmaxt  at  each  tr.  When  >  2,  Ci^maxi  and  Ci^rma  be¬ 
come  independent  of  <r,  but  the  magnitudes  of  Ci^maxi^ 
Cm, maxi  and  CM^rms  for  =  4  are  smaller  than  those 
for  cr  =  2  and  3.  When  <r  approaches  zero,  all  the  coeffi¬ 
cients  tend  to  be  proportional  to  ((rvmaxt)  or  {<rvmax) 
which  is  proportional  to  the  circulation  of  the  vortex 
tube.  For  example,  CL,rms  la  expressed  by 

Cl  ,rmj  ^1  ^  ^  ^  ^  4 

=  C2  Vmaxt  0.25  <  ^  <  2, 

0.75>n>0;3,  (17) 

where  the  constant  ct  =  1  and  C2  —  0.7,  and  n  de¬ 
pends  on  <r  and  should  approach  unity  as  cr  reaches 
zero.  CL^maxl  >  CLjmax7%  Cm  ,ma«}  and  CM,rnu  are  also 
expressed  by  equation  (17)  with  ci  =  1.1  and  =  0,88, 

Cl  =  ~2  and  cj  =  —1.65,  ci  =  0.13  and  ca  =  0.11,  and 
Cl  =  0.053  and  cj  =  0.04,  respectively.  The  time  aver¬ 
aged  value  of  the  deviation  of  the  drag  coefficient  from 
that  of  the  axisymmetric  flow  past  a  sphere  for  all  val¬ 
ues  of  (r  is  nearly  zero  (©(lO^'*)). 

Note  that  CL,max2i  Cm, maxi  and  CM,rm§  for  cr  =  4 
are,  respectively,  smaJler  than  those  for  a  =  2  amd  3, 
and  the  reason  is  explauned  as  follows.  When  the  initial 
size  of  the  vortex  core  is  considerably  larger  than  the 
sphere  size  (<r  >  4),  the  effect  of  the  shear  flow  (induced 
by  the  passage  of  the  vortex  tube)  across  the  sphere 
diminishes.  As  a  result,  the  magnitude  of  CL,max2  de¬ 
creases  and  approaches  Ci^maxi  as  (r  »  1.  We  ex¬ 
plained  in  section  3.2.1  that  the  magnitude  of  CL,max2 
depends  on  the  combined  effect  of  the  downwash  and 
the  shear  flow  across  the  sphere  due  to  the  vortex  tube. 
Cm, max  and  CM,rma  also  decreases  for  the  same  reason. 

In  addition,  when  the  initial  size  of  the  vortex  core  is 
larger  than  the  sphere  size,  the  effect  of  the  wake  be¬ 
hind  the  sphere  on  the  vortex  tube  diminishes.  As  a 
consequence,  the  magnitudes  of  the  lift  and  moment 


coefficients  decay  slowly  towards  zero  after  they  peak 
near  the  time  of  passage  of  the  pair  of  vortex  tubes  by 
the  sphere. 

In  sununary,  the  maximum  positive  lift  coefficient 
and  the  rms  lift  coefficient  depend  only  on  circulation 
at  small  values  of  <t  while  they  depend  only  on  Vmaxt 
(and  not  (t)  at  large  values  of  a.  For  mid-range  values 
of  cr,  they  depend  on  both  (t  and  Vmaxt  (or  equivalently 
both  <T  and  F/). 

Now,  the  effects  of  the  offset  distance  on  the  flow 
field  are  investigated  by  varying  dojj  for  Re  =  100  and 

=  4.  The  computation  was  performed  for  /  =:  ±1, 
±2,  ±3,  ±4,  and  0  in  addition  to  the  base  case  d®//  = 
±1.5.  Note  that  the  case  of  do/ /  =  0  corresponds  to  the 
interaction  between  a  single  vortex  tube  and  a  sphere. 

It  is  found  that  C L,maxXt  Cl  ,max7t  Cl  ,rm$i  Cm  ,maxt 
and  CMrmt  for  «ach  do//  are  linearly  proportional  to 
Vmaxt  as  in  the  case  of  do//  =  ±1.5.  The  triangular 
symbols  in  Figure  6  show  CL,rm*  as  a  function  of  |do//| 
for  Re  =  100  amd  =  4  while  the  maximum  fluctuar 
tion  velocity  (or  the  circulation)  of  each  vortex  tube  is 
kept  as  a  constant,  Vmax  =  0.2.  The  figure  shows  that 
CL,rm»  decays  rapidly  as  ldo//|  >  0  for  Vmax  =  0.2. 
On  the  other  hand,  the  circular  symbols  in  Figure  6 
show  CL,rmt  as  a  function  of  |do//l  for  Re  =:  100  and 

=  4  while  the  total  maximum  fluctuation  velocity  in¬ 
duced  by  the  two  vortex  tubes  is  kept  as  a  constant, 
Vmaxt  =  0.2.  They  also  show  that  the  magnitudes  of 
the  rms  lift  coefficients  for  d^//  =  ±1,  ±1.5,  ±2,  ±3, 
and  ±4  are  close  to  that  for  do//  =  0.  After  exam¬ 
ining  the  effect  of  offSset  distance  for  tr  =  1  and  2,  we 
found  that  the  equation  of  the  rms  lift  coefficient  for  a 
single  vortex  tube  can  be  applicable  to  that  for  a  pair 
of  vortex  tubes  when  the  separation  distance  between 
their  centers  is  less  than  2  ^/^  vortex  tube  diameter 
for  Re  =  100  if  Vmaxt  is  used  instead  of  Vmax  in  the 
equation  of  the  rms  lift  coefficient.  The  behaviors  of 
Cl, maxi  and  CL,max2  as  a  function  of  jdo//|  are  similar 

to  that  of  Cl, rm#. 

The  magnitude  of  the  rms  moment  coefficient  as  a 
function  of  |do//|  decays  faster  than  that  of  the  rms 
lift  coefficient  as  d<,//  >  0,  and  the  equation  of  the 
rms  moment  coefficient  for  a  single  vortex  tube  can  be 
applicable  to  that  for  a  pair  of  vortex  tubes  when  the 
separation  distance  between  their  centers  is  less  than 
^/a  vortex  tube  diameter  for  Re  =  100  if  Vmaxt  is  used 
instead  of  Vmax  in  the  equation  of  the  rms  moment 
coefficient.  The  behavior  of  Cm, max  ^  ^  function  of 
do//  is  similar  to  that  of  CM,rm$- 

3.2.3  Effects  of  Reynolds  number 

Similar  computations  to  those  in  section  3.2.2  are 
made  for  four  different  Reynolds  numbers  in  the  range 
of  20  <  i2e  <  80,  do//  =  ±1.5,  and  (T  =  1  in  addition 
to  the  base  case  Re  =  100. 

CL,maxliCL,max2iCL,rmaiC M ,maxi  ^d  CM,rma 
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linearly  proportional  to  Vmaxt  for  each  Reynolds  num¬ 
ber  as  in  the  case  of  =  100.  Figure  7  shows 
and  Ci^rmt  as  a  function  of  Reynolds  num¬ 
ber  for  dojj  =  ±1.5  and  <r  =  1  with  Vmaxt  =  0,295. 
Cx,m<jxi  increases  with  a  negative  power  of  Reynolds 
number  as  Reynolds  number  decreases  and  is  expressed 
by 

Cx,maxl  =  A  Vrtxaxi  ,  (18) 

where  the  constant  A  =  4.6  and  P  =  —0.37. 

(^L,max2  is  not  Sensitive  to  the  change  of  the  Reynolds 
numbers  and  slowly  increase  as  Reynolds  number  de¬ 
creases.  When  Reynolds  number  decreases,  the  effect 
of  wake  behind  the  sphere  decreases,  and  thus  the  lift 
coefficient  decays  slowly  to  zero  after  it  attains  the  max¬ 
imum  negative  value.  Ci^rmt  increases  with  a  negative 
power  of  Reynolds  number  as  Reynolds  number  de¬ 
creases  and  is  expressed  by  equation  (18)  with  A  =  3.5 
and  P  =  -0.35.  The  maximum  moment  coefficient 
and  the  rms  moment  coefficient  are  affected  by  only 
the  viscous  effect  and  increases  with  a  negative  power 
of  Reynolds  number  as  Reynolds  number  decreases. 
Csf  ,max  Cm, rms  follow  the  form  of  equation  (18) 
with  A  =  3.26  and  P  =  —0.75  for  the  former,  and 
A  =  1.43  and  P  =  —0.77  for  the  latter  for  <r  =  1. 

The  effect  of  Reynolds  number  (20  <  iEe  <  80)  on  the 
flow  field  is  also  investigated  for  the  vortex  size  larger 

than  O'  =  1  (2  <  (T  <  4)  and  do//  =  ±L5.  Cc^maxi 

and  Ci^rms  are  linearly  proportional  only  to  Vmaxt  and 
independent  of  <r  when  (T  >  2  at  fixed  Reynolds  num¬ 
ber  as  in  the  case  of  Re  =  100.  Cl, maxi  and  Cl, rms 
increase  with  a  negative  power  of  Reynolds  number  as 
Reynolds  number  decreases  and  follow  the  form  of  equa¬ 
tion  (18)  with  ^  =  8.9  and  P  =  -0.45  for  the  for¬ 
mer,  and  =  8.1  and  P  =  —0,45  for  the  latter  for 
2  ^  £  4.  Cm, max  and  Cm, rms  also  increase  with  a 

negative  power  of  Reynolds  number  as  Reynolds  num¬ 
ber  decreases  and  follow  the  form  of  equation  (18)  with 
A  =  5.5  and  P  =  —0.83  for  the  former,  and  i4  =  3.1 
and  P  =  -0,88  for  the  latter  for  2  <  <r  <  4. 

Now,  the  effect  of  the  offset  distance  is  determined  for 
20  <  Pe  <  80,  in  addition  to  the  base  case  Re  =  100. 

The  triangular  symbols  in  Figure  8  show  Cl, rms  as 
a  function  of  jd<,//)  for  Re  =  20  and  (t  =  4  while  the 
maximum  fluctuation  velocity  (or  the  circulation)  of 
each  vortex  tube  is  kept  as  a  constant,  Vmax  =  0.2.  The 
figure  shows  that  Cl, rms  decays  fast  as  |do//|  >  0  for 
Vmax  =  0.2.  On  the  other  hand,  the  circular  symbols  in 
Figure  8  show  Cl, rms  as  a  function  of  \dof  / 1  for  Pe  =  20 
and  (T  =  4  while  the  total  maximum  fluctuation  velocity 
induced  by  the  two  vortex  tubes  is  kept  as  a  constant, 
Vmaxt  —  0.2.  They  also  show  that  the  magnitudes  of  the 
rms  lift  coefficients  for  d^jj  =  ±2  and  ±4  are  close  to 
that  for  d<,//  =  0.  The  results  for  the  range  of  (x  values 
indicate  that  the  equation  of  the  rms  lift  coefficient  for 
a  single  vortex  tube  can  be  applicable  to  that  for  a  pair 
of  vortex  tubes  when  the  separation  distance  between 
bheir  centers  is  less  than  2  vortex  tube  diameter  for 


Pc  =  20  if  Vmaxt  is  used  instead  of  Vmax  in  the  equation  ' 
of  the  rms  lift  coefficient.  The  same  result  as  above  was 
obtained  at  different  Reynolds  numbers,  Re  =:  40,  60, 

Md  80.  The  behaviors  of  CL,maxi  and  Cl^ox^  as  a 
function  of  Ido//|  are  similar  to  that  of 

The  magnitude  of  the  rms  moment  co^cient  as  a 
function  of  \dojj\  decays  faster  than  that  of  the  rms  lift 
coefficient  as  d^fj  becomes  large,  and  the  rms  moment 
coefficient  for  a  single  vortex  tube  can  be  applicable 
to  that  for  a  pair  of  vortex  tubes  when  the  separation 
distance  between  their  centers  is  less  than  vortex 
tube  diameter  if  Vmaxt  is  used  instead  of  Vmax  in  the 
equation  of  the  rms  moment  coefficient.  The  behavior 
of  Cm, max  as  a  function  of  |do//|  is  similar  to  that  of 

Cm, rms  • 

In  simmary,  comparison  of  the  results  from  this  sec¬ 
tion  with  those  from  the  previous  section  shows  that 
the  range  of  the  offset  distance  for  which  the  equa¬ 
tions  for  a  single  vortex  tube  can  be  applicable  to  those 
for  a  pair  of  vortex  tubes  does  not  change  as  a  func¬ 
tion  of  Reynolds  number  for  20  <  Pe  <  100,  Cl  rms 

Cm, rms  for  a  single  vortex  tube  which  were  ob¬ 
tained  by  Kim  ei  al  are  expressed  by  equation  (18)  with 
A  =  8.1  and  P  =  —0.45  for  the  former,  and  A  =  3.1 
and  P  =  — 0.88.  The  rms  lift  coefficient  is  written  again 
here  for  later  use. 

CL,rms  =  8.1  Urnart  ,  2  <  <7  <  4  (19) 

3.3  Interactions  of  a  pair  of  vortex  tubes  of  opposite 
rotation  and  a  sphere 

We  consider  the  same  initial  flow  geometry  and  pa¬ 
rameters  as  those  in  section  3.2  but  for  a  pair  of  vortex 
tubes  of  opposite  rotation.  The  base  case  calculation 
is  that  of  Pe  =  100,  ±1.5,  and  <r  =  1,  Note 

that  the  lift  and  torque  on  the  sphere  are  zero  due  to 
the  flow  symmetry  in  upper  and  lower  regions  of  the 
sphere 


Figure  9  shows  the  drag  coefficients  of  the  sphere 
as  a  function  of  time  for  Pe  =  100,  d^//  =  ±1.5,  and 
(T-l.  The  drag  coefficients  are  obtain^  with  four  dif¬ 
ferent  total  maximum  fluctuation  velocities  due  to  the 
vortex  tubes,  Vmaxt  =  0.185,  0.369,  0.554,  and  0.738 
(Vm4ix  =  0.1,  0.2,  0.3,  and  0.4).  The  temporal  behav¬ 
ior  of  the  drag  coefficients  is  different  from  that  of  the 
case  of  the  pair  of  vortex  tubes  of  like  rotation.  The 
time-averaged  value  of  the  deviation  of  the  drag  coeffi¬ 
cient  from  that  of  the  axisymmetric  flow  past  a  sphere 
for  all  values  of  Vmaxt  is  not  negligible  and  increased 
linearly  as  Vmaxt  increases.  When  the  top  and  bottom 
vortex  tubes  have  positive  and  negative  circulations,  re¬ 
spectively,  the  induced  velocity  due  to  the  vortex  tubes 
adds  its  magnitude  to  the  base  flow  along  the  stagna¬ 
tion  streamline,  and  so  the  dynamic  pressure  ahead  of 


3*3.1  Tubes  of  top-positive  and 
bottom-negative  circulations 
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the  front  stagnation  point  becomes  higher  than  that  of 
the  axisymmetric  flow  past  a  sphere.  This  causes  the 
pressure  at  the  stagnation  point  and  the  shear  stresses 
in  the  upper  and  lower  left  regions  to  be  higher  than 
those  of  the  axisymmetric  flow  past  a  sphere.  As  a 
consequence,  the  drag  is  increased. 

3,3«2  Tubes  of  top-negative  and 
bottom-positive  circulations 

Figure  10  shows  the  drag  coefficients  of  the  sphere  as 
a  function  of  time  for  the  same  parameters  as  used  in 
section  3.3.1.  The  drag  coefficients  are  obtained  with 
four  different  total  maximum  fluctuation  velocities  due 
to  the  vortex  tubes,  Vmaxt  =  0.185,  0.369,  0.554,  and 
0.738  {vmas  =  0-1,  0.2,  0.3,  and  0.4).  The  temporal 
behavior  of  the  drag  coefficients  is  different  from  that 
of  the  case  of  the  pair  of  vortex  tubes  of  like  rotation. 
The  time  averaged  value  of  the  deviation  of  the  drag 
coefficient  from  that  of  the  axisymmetric  flow  past  a 
sphere  for  all  values  of  Vmaxt  in  not  negligible  and  de¬ 
creased  linearly  as  Vmaxt  increases.  When  the  top  and 
bottom  vortex  tubes  have  negative  and  positive  circular 
tions,  respectively,  the  induced  velocity  due  to  the  vor¬ 
tex  tubes  subtracts  its  magnitude  from  the  base  flow 
adong  the  stagnation  stre2unline,  and  so  the  dynamic 
pressure  ahead  of  the  front  stagnation  point  becomes 
lower  than  that  of  the  axisymmetric  flow  past  a  sphere. 
This  causes  the  pressure  at  the  stagnation  point  and 
the  shear  stresses  in  the  upper  and  lower  left  regions 
to  be  lower  than  those  of  the  axisymmetric  flpw  past  a 
sphere.  As  a  consequence,  the  drag  is  decreased. 

4.  Conclusions 

In  order  to  understand  better  the  physics  of  interac¬ 
tion  between  a  particle  and  the  turbulent  carrier  flow, 
the  unsteady,  three-dimensional,  incompressible,  vis¬ 
cous  flow  interactions  between  a  pair  of  vortex  tubes 
advected  by  a  uniform  free  stream  and  a  spherical  par¬ 
ticle  suddenly  placed  and  held  fixed  in  space  is  investi¬ 
gated  numerically  for  a  range  of  particle  Reynolds  num¬ 
ber  20  <  iZe  <  100.  The  counter-clockwise  rotating  a 
pair  of  cylindrical  vortex  tubes  are  initially  located  ten 
radii  upstream  from  the  center  of  the  sphere. 

A  summery  of  the  findings  and  their  applications  is 
provided  as  follow. 

(i)  The  effects  of  the  size  and  the  offset  distance  of 
the  pair  of  vortex  tubes  on  the  flow  field  are  examined 
for  20  <  Re  <  100.  The  lift  and  moment  coefficients 
are  found  to  be  linearly  proportional  to  the  maximum 
fluctuation  velocity  (v^oxt)  induced  by  the  pair  of  vor¬ 
tex  tubes  at  given  size  (a)  and  offset  distance  of  the 
vortex  tube,  and  the  rms  lift  coefficient  depends  only 
on  Vmor#  but  independent  of  <t  when  >  2.  Further¬ 
more,  the  equations  for  the  lift  and  moment  coefficients 
of  the  sphere  for  a  single  vortex  tube  are  applicable  to 
those  for  a  pair  of  vortex  tubes  when  the  separation 


distance  between  their  centers  is  less  than  2^/?  vor¬ 
tex  tube  diameter  for  the  lift  coefficients  and  less  than 
y/r  vortex  tube  diameter  for  the  moment  coefficients  if 
Vmaxt  is  used  instead  of  Vmax  in  the  equations.  These 
separation  distances  do  not  vary  for  a  range  of  Reynolds 
number  20  <  Re  <  100. 

(ii)  The  results  in  (i)  can  be  applied  to  turbulent  flows 
containing  small  concentration  of  particles  in  order  to 
obtain  the  rms  lift  force  on  a  peirticle.  A  turbulent  flow 
possesses  a  wide  spectrum  of  eddy  sizes.  The  large  ed¬ 
dies  contain  moat  of  the  turbulent  kinetic  energy  and 
produce  high  velocity  fluctuations,  and  so  they  are  re¬ 
sponsible  for  the  dispersion  of  particles.  The  particle 
size,  at  the  extremes,  may  be  comparable  to  either  the 
integral  length  scale  or  to  the  Kolmogorov  length  scale.. 
When  the  size  of  particle  is  comparable  to  the  integral 
length  scale,  the  rms  lift  coefficient  of  the  particle  is  ob¬ 
tained  by  equation  (19).  Furthermore,  the  results  tend 
to  support  the  idea  that  equation  (19)  would  be  appli¬ 
cable  to  the  case  of  an  eddy  being  much  larger  than  the 
particle.  Thus,  when  the  size  of  particle  is  comparable 
to  the  Kolmogorov  length  scale,  the  rms  lift  coefficient 
of  the  particle  can  be  obtained  approximately  by  equar 
tion  (19),  where  Vmaxt  is  the  maximum  velocity  fluctu¬ 
ation  induced  by  eddies  with  the  integral  length  scale. 
The  time  during  which  the  particle  is  influenced  by  the 
eddy  is  of  the  order  of  the  eddy  lift  time. 

The  deflection  of  the  particle  path  will  depend  on 
the  magnitude  of  the  rms  lift  coefficient  cind  the  ratio, 
Pr,  of  the  particle  density  to  that  of  the  carrier  fluid 
(Cx  =  \prA,  where  A  is  the  dimensionless  accelera¬ 
tion  of  the  particle  due  to  the  lift  force).  This  result 
provides  a  simple  method  to  estimate  the  deflection  of 
particle  trajectory  in  the  dilute  particle-laden  turbu¬ 
lent  flow.  Equation  (19)  and  the  nondimensionalized 
Newton’s  second  law  show  that  the  deflection  increases 
slowly  as  Reynolds  number  decreases. 

(iii)  The  magnitude  of  the  rms  moment  coefficient  of 
the  particle  is  one  order  of  magnitude  less  than  that 
of  the  rms  lift  coefficient  when  Re  >  20.  Furthermore, 
when  the  initial  size  of  the  vortex  core  is  considerably 
larger  than  the  sphere  size  (o’  >  4),  the  effect  of  the 
shear  flow  (induced  by  the  passage  of  the  vortex  tube) 
aw:ross  the  sphere  diminishes  and  the  torque  on  the  par¬ 
ticle  decreases.  Thus,  the  torque  on  the  particle  might 
be  negligible  in  many  applications. 

(iv)  When  the  top  and  bottom  vortex  tubes  have 
positive  and  negative  circulations,  respectively,  the  in¬ 
duced  velocity  due  to  the  vortex  tubes  adds  its  mag¬ 
nitude  to  the  base  flow  along  the  stagnation  stream¬ 
line.  This  causes  the  pressure  at  the  stagnation  point 
and  the  she2ur  stresses  in  the  upper  and  lower  left  re¬ 
gions  to  be  higher  than  those  of  the  axisymmetric  flow 
past  a  sphere.  As  a  consequence,  the  drag  is  increased. 
On  the  other  hand,  When  the  top  and  bottom  vortex 
tubes  have  negative  and  positive  circulations,  respec¬ 
tively,  the  induced  velocity  due  to  the  vortex  tubes 
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subtracts  its  magnitude  from  the  base  flow  along  the 
stagnation  streamline.  This  causes  the  pressure  at  the 
stagnation  point  and  the  shear  stresses  in  the  upper 
and  lower  left  regions  to  be  lower  than  those  of  the  ax- 
isymmetric  flow  past  a  sphere.  As  a  consequence,  the 
drag  is  decreased. 
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Figure  1.  Flow  geometry  and  coordinates 


Figure  2.  Contour  lines  of  y-component  vorticity 
in  the  principal  plane  at  (a)  f  =r  1,  (b)  6, 
(c)  10,  (d)  15,  (e)  21,  and  (f)  30  for 
Re  =:  100,  doff  =  0  and  <r  =:  1 
with  Vmest  =  0.59. 
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Figure  3.  Lift  coefficients  of  the  sphere  as  a 
function  of  time  and  Vmaxt  ^or 
Rc  =  100,  doff  =  ±1.5,  (T  =  1. 

*  denotes  the  case  of  a  single  vortex. 
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Figure  5.  Drag  coefficients  of  the  sphere  under 
the  same  conditions  of  figure  3. 
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Figure  6.  Rms  lift  coefficients  of  the  sphere  as  a 

function  of  \dof  / 1  for  iZe  =  100  and  fr  =  4. 


Figure  4.  Moment  coefficients  of  the  sphere  under 
the  same  conditions  of  figure  3. 

*  denotes  the  case  of  a  single  vortex. 
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Figure  7.  Maximum  positive  lift  coefBcient  and 
rms  lift  coefficient  as  a  function  of 
Reynolds  number  for  d^ff  =  ±1.5 
and  cr  =  1  with  Vmaxt  =  0.295. 


Figure  9.  Drag  coefficients  of  the  sphere  as  a 
function  of  time  and  Vmaxt  for 
Re  =  100,  doff  =  ±1.5,  and  0"=  1 
with  top-positive  and  bottom-negative 
circulations. 


Figure  8.  Rms  lift  coefficients  of  the  sphere  as  a 

u.a  -  '  ' 

\  ! 

function  of  \do/f\  for  Re  =  20  and  cr  =  4. 

\  / 

\  / 

0  5  10  15  20 

t 

Figure  10.  Drag  coefficients  of  the  sphere  as  a 

Bi-C'" 

function  of  time  and  Vmast  for 

pur:' 

Re  =  100,  doff  =  ±1.5,  amd  (T  =  1 
with  top-negative  and  bottom-positive 
circulations. 

